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Executive  Summary 

The  objective  of  our  SEED  project  was  to  derive  a  method  to  trade  habitat  that  incorporates  the 
role  of  landscape  connectivity  for  population  viability,  and  level  of  scientific  uncertainty 
regarding  landscape-scale  processes.  Habitat  trading  provides  an  opportunity  to  change  the 
spatial  allocation  of  protected  habitat  to  relieve  training  restrictions,  but  may  also  threaten  the 
long-term  viability  of  some  populations.  Some  species  are  able  to  achieve  regional  persistence 
through  the  exchange  of  individuals  among  local  breeding  groups  while  maintaining  a  certain 
level  of  independence  in  productivity  at  local  scales,  thereby  avoiding  spatial  correlation  in 
demographic  cycles  that  could  lead  to  regional  extinction  (Hanski  and  Gaggiotti  2004). 

Therefore,  landscape  connectivity,  or  the  extent  to  which  landscape  patterns  facilitate  movement 
of  individuals,  may  contribute  substantially  to  population  viability  for  many  species  found  on 
DoD  installations. 

Our  study  focused  on  the  Red-Cockaded  Woodpecker  (RCW)  because  its  presence  on  many 
DoD  facilities  in  the  southeastern  US  has  led  to  training  restrictions.  Further,  there  is  a  scientific 
consensus  that  habitat  connectivity  is  critical  for  population  persistence  (Cox  et  al.  2001;  Schiegg 
et  al.  2002).  RCWs  are  well  studied  due  to  their  cooperative  mating  system  and  tendency  to  pass 
breeding  territories  to  the  next  generation.  A  breeding  pair  inhabits  a  territory  and  often  shares 
the  habitat  with  adult  helpers  who  contribute  to  the  defense  of  the  territory  and  recruitment  of 
offspring  (Conner  et  al.  2001). 

We  chose  to  work  on  U.S.  Marine  Corps  Base  Camp  Lejeune  (MCBCL)  because  the  natural 
resource  managers  expressed  an  interest  in  helping  to  restore  RCW  habitat  off-base  to  relieve 
training  restrictions  and  prevent  urban  encroachment.  MCBCL  is  home  to  the  Marine’s  2nd 
Expeditionary  Force  and  is  located  on  the  southern  coast  of  North  Carolina  in  Onslow  County. 
MCBCL  has  already  completed  cooperative  purchases  of  five  off-site  parcels  through 
partnerships  with  North  Carolina  Wildlife  Resources  Commission  (NCWRC)  and  North  Carolina 
Coastal  Land  Trust  (NCCLT),  in  accordance  with  the  U.S.  Army’s  Compatible  Use  Buffer 
Program  and  recommendations  of  Onslow  County’s  Joint  Landuse  Study.  The  cooperative 
purchase  of  off-base  properties  to  prevent  residential  or  commercial  development  is  often 
referred  to  as  Encroachment  Partnering. 

Determining  how  Encroachment  Partnering  will  lead  to  regulatory  relief  on-base  is  an  important 
area  of  research  because  habitat  trading  programs  will  change  how  habitats  are  allocated  for 
small  and  sometimes  declining  populations.  Uncertainty  regarding  critical  biological  processes 
required  to  maintain  persistent  subdivided  populations  should  be  the  focus  of  any  scientifically 
defensible  habitat  trading  program.  We  addressed  two  critical  uncertainties.  First,  we  were 
uncertain  how  differences  in  vegetative  composition  of  territories  affect  reproductive  success 
(James  et  al.  2001;  Walters  et  al.  2002).  Second,  we  were  uncertain  how  land  use  patterns 
between  RCW  territories  affect  the  probability  of  successful  dispersal  (J.  Walters,  personal 
communication). 

This  project  derived  an  approach  to  integrate  uncertainty  analysis  of  critical  biological  processes 
with  a  landscape-scale  tradable  credit  system  designed  to  minimize  habitat  fragmentation  effects 
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on  population  viability.  Our  project  focuses  on  the  development  and  application  of  spatially- 
explicit  population  models  (SEPMs)  to  inform  habitat  trading  decisions.  SEPMs  represent  a 
hypothesis  regarding  how  a  population  may  respond  to  changing  land  use  patterns.  Our  project 
can  be  divided  into  three  phases.  The  first  was  to  use  Pattern-Oriented  Modeling  to  calibrate  the 
ROW  SEPM  for  MCBCL  while  testing  a  new  submodel  to  describe  dispersal  behaviors.  The 
second  was  to  use  all  SEPMs  capable  of  reproducing  multiple  patterns  to  evaluate  habitat  trades 
on  and  near  MCBCL  using  the  generalized  landscape-scale  accounting  system  Landscape 
Equivalency  Analysis.  The  third  phase  was  to  develop  a  new  recruitment  submodel  that  may  be 
used  to  extend  the  RCW  SEPM  to  areas  that  lack  historical  demographic  data  available  on 
MCBCL  but  would  still  like  to  participate  in  habitat  trading  with  DoD. 

Phase  1:  Uncertainty  Regarding  Dispersal 

Despite  a  tremendous  amount  of  knowledge  regarding  the  mating  system  of  RCWs,  we  are  still 
uncertain  how  RCWs  move  through  a  heterogeneous  landscape.  We  tested  a  new  dispersal 
submodel  using  Pattern-Oriented  Modeling  (POM).  POM  is  a  structured  approach  for 
understanding  the  underlying  processes  that  generate  the  patterns  we  observe  in  nature  (Grimm 
et  al.  2005).  POM  tests  the  ability  of  several  thousand  SEPMs,  each  with  a  different 
parameterization,  to  reproduce  patterns  observed  in  nature.  Often  just  a  few  SEPMs  satisfy  all 
pattern  fulfillment  criteria.  POM  provides  an  objective  and  quantitative  approach  for  including 
multiple  models  and  consequently  uncertainty  into  model-based  environmental  decision  making. 

We  compiled  pre-existing  data  on  land  cover  patterns  and  patterns  describing  the  RCW 
population  on  MCBCL.  Next,  the  results  of  100,000  different  SEPMs,  each  describing  different 
dispersal  behaviors,  were  compared  with  RCW  patterns  observed  on  MCBCL.  Only  8  models 
were  capable  of  reproducing  all  the  patterns — a  significant  reduction  in  uncertainty.  For 
parameters  relating  to  floater  survival  and  distance  at  which  a  breeding  vacancy  can  be  detected, 
we  were  able  to  reduce  parameter  uncertainty  by  54-92%.  However,  for  parameters  relating  to 
the  attractiveness  of  different  land  cover  types  for  dispersal,  or  conversely  landscape  friction,  we 
were  unable  to  reduce  parameter  uncertainty.  This  may  have  resulted  because  existing  territory 
locations  are  highly  clumped  on  MCBCL  and  floaters  do  not  have  to  interact  with  much  of  the 
landscape  in  order  to  find  a  breeding  vacancy.  POM  provides  a  method  of  improving  SEPMs 
without  spending  more  money  on  field  data  collection  and  can  be  used  to  prioritize  which 
patterns  in  nature  are  most  useful  for  reducing  prediction  uncertainty  (Wiegand  et  al.  2004). 

Phase  2:  Habitat  Trading 

Habitat  patches  often  differ  in  conservation  value  due  to  their  size,  shape,  vegetative 
composition,  and  location  relative  to  other  land  cover  types  in  the  landscape.  We  devised  the 
Landscape  Equivalency  Analysis  (LEA)  to  determine  if  habitat  patches  traded  make  equivalent 
contributions  to  rates  of  recruitment  within  breeding  groups  and  rates  of  dispersal  among 
breeding  groups  (Bruggeman  et  al.  2005).  LEA  assigns  conservation  value  to  habitat  at  the  local 
scale  based  on  its  contribution  to  abundance  and  genetic  variance  estimated  at  the  regional  scale. 

The  uniqueness  of  this  approach  lies  in  the  ability  of  LEA  to  incorporate  the  unequal  contribution  of 
habitat  patches  and  areas  between  patches  into  the  estimate  of  conservation  value.  A  tradable  credit 
system  that  accounts  for  differences  in  conservation  value  among  habitat  patches  is  more  likely 
to  achieve  sustainable  populations  at  the  least  cost. 
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In  Phase  2  we  used  Decision  Analysis  which  provides  a  structured  approach  for  comparing 
alternative  choices  in  the  face  of  scientific  uncertainty.  Alternative  mitigation  choices  were 
evaluated  using  the  eight  SEPMs  derived  from  POM  and  the  Letcher  SEPM,  the  model  currently 
used  by  DoD  and  DOI  to  inform  ROW  management  decisions.  For  all  habitat  trades  and  models 
considered,  restoring  ROW  habitat  outside  the  Main  Base  to  compensate  for  the  loss  of  clusters 
inside  the  Main  Base  increased  habitat  fragmentation  effects.  The  Letcher  model  predicted  the 
smallest  increase  in  fragmentation  effects.  Valuing  credits  based  on  habitat  occupancy,  or  the 
number  of  Potential  Breeding  Groups  (PBGs),  indicated  that  restoring  habitat  on  all  parcels 
outside  of  the  Main  Base,  including  Lanier,  Beck,  Allen  and  the  Greater  Sandy  Run  Area,  is  the 
best  decision.  In  contrast,  valuing  credits  based  on  changes  in  rates  of  recruitment  within 
clusters  and  rates  of  migration  among  clusters  indicated  that  restoring  only  the  Lanier  tract  is  the 
best  habitat  trade.  The  Lanier  tract  was  the  smallest  parcel  considered,  but  provided  the  greatest 
connectivity  with  RCW  clusters  on  the  Main  Base.  Decision  analysis  indicated  that  despite  the 
remaining  uncertainty  in  dispersal  behaviors,  a  further  reduction  in  uncertainty  will  likely  not 
change  the  decision  to  restore  only  the  Lanier  tract  instead  of  all  possible  locations. 

LEA  estimates  the  financial  value  of  tradable  credits  based  on  the  in-kind  replacement  costs  of 
ecological  services  that  differ  in  sensitivity  to  changes  in  habitat  area  and  connectivity. 

Therefore,  LEA  provides  a  nexus  for  ecological-economic  analysis  of  management  decisions  at  a 
landscape  scale.  Even  though  no  trade  considered  prevented  habitat  fragmentation  effects,  each 
trade  did  increase  habitat  occupancy.  For  discussion  purposes  we  calculated  the  value  of  a  LEA 
credit  for  Potential  Breeding  Groups  and  discussed  how  DoD  could  detennine  when  it  is  more 
cost-effective  to  forego  a  trade  to  collect  more  data  and  update  model  assumptions  to  re-evaluate 
trading  choices. 

Phase  3:  Habitat-Fitness  Relationships 

In  an  effort  to  increase  the  spatial  scale  of  our  model,  we  investigated  statistical  relationships 
between  pre-existing  data  describing  foraging  habitat  on  MCBCL  and  the  number  of  fledglings 
produced  per  cluster.  Early  in  the  project  we  detennined  that  data  describing  RCW  demography 
required  by  the  current  SEPM  is  not  available  across  the  Onslow  Bight  Landscape.  The  SEPM 
has  been  constructed  to  work  with  the  RCW  monitoring  programs  implemented  on  MCBCL  and 
Fort  Bragg,  which  have  contributed  much  to  management  objectives  at  smaller  scales  and  the 
field  of  landscape  ecology.  We  argue  that  DoD’s  goal  of  using  Encroachment  Partnering  to 
contribute  to  regulatory  relief  requires  updating  model  structure  to  work  at  broader  spatial  scales. 

The  current  model  estimates  the  number  of  fledglings  per  group  (group  fitness)  based  on  the  ages 
of  parents  and  number  of  helpers.  While  this  approach  is  based  on  a  statistically  significant 
regression  equation  (p<0.0001,  R"  =  0.08),  it  requires  knowing  the  age  structure  of  the  breeding 
population,  which  may  be  determined  using  long-term  bird  banding  studies.  Replicating  this 
effort  across  the  Onslow  Bight  Landscape  would  be  difficult. 

Previous  researchers  have  attempted  to  relate  group  fitness  to  vegetative  characteristics  near 
cavity  trees  (James  et  al.  2001;  Convery  2002;  Walters  et  al.  2002).  While  these  efforts  have 
observed  that  the  number  of  adults  in  a  group  is  related  to  vegetative  characteristics,  statistically 
significant  relationships  with  fitness  were  often  not  observed.  However,  Convery  (2002)  did 
successfully  relate  group  fitness  to  habitat  attributes  in  23  territories  on  MCBCL.  Convery 


3 


Final  Draft  v2 


Bruggeman,  Jones,  and  Wiegand 


(2002)  supplemented  pre-existing  data  on  forest  structure  with  surveys  of  ground  cover  within  a 
group’s  home  range.  Defining  a  group’s  home  range  and  characterizing  ground  cover  is  labor 
intensive  and  these  efforts  have  not  yet  been  replicated  for  other  groups  on  MCBCL. 

In  this  phase,  we  explored  habitat-fitness  relationships  using  a  moving  window  approach  and 
incorporating  yearly  variation  in  recruitment.  The  moving  window  approach  was  used  because 
we  will  often  lack  the  resources  to  define  home  ranges.  Instead  we  characterized  the  forest 
structure  within  500  m,  1,000  m,  and  3,500  m  of  the  cavity  tress.  The  RCW  Foraging  Habitat 
Survey  conducted  on  MCBCL  from  1998  to  2000  was  used  to  estimate  the  density  of  pines  and 
hardwoods  within  these  windows.  Further,  field  observations  suggested  that  variation  in  rainfall 
may  affect  recruitment.  To  conduct  an  initial  screen  of  this  hypothesis,  we  included  a  year  effect 
to  determine  if  year  to  year  variation  in  environmental  factors  may  be  an  explanatory  variable  for 
fitness.  Selection  among  alternative  models  was  performed  using  the  Akaike  Information 
Criterion  (AIC)  to  determine  the  most  parsimonious  model  capable  of  predicting  recruitment 
(Burnham  and  Anderson  2002). 

We  were  unable  to  derive  a  statistically  significant  regression  predicting  group  fitness  when  only 
habitat  variables  were  included  (p=0.063,  R2=0.03),  but  these  results  did  suggest  that  forest 
structure  at  the  1,000  m  and  3,500  m  might  affect  recruitment.  Next,  we  included  the  influence 
of  helpers,  which  improved  model  fit  (p=0.0014,  R2=0.09).  Variables  describing  forest  structure 
within  3,500  m  were  still  retained  in  this  regression  equation.  Inclusion  of  a  year  effect  further 
improved  model  fit  (p=0.0001,  R2=0.12),  and  increased  the  significance  of  the  habitat  variables. 
We  found  that  the  year  effect  was  more  important  than  habitat  features  but  less  important  than 
the  number  of  helpers.  Researchers  have  argued  that  models  with  AIC  scores  within  2  units 
should  be  treated  as  equivalent  (Richards  2005).  Therefore,  we  further  simplified  the  models 
until  they  were  within  2  AIC  units  of  the  model  with  the  lowest  AIC  score.  This  approach 
indicated  that  the  number  of  fledglings  could  be  predicted  based  on  the  number  of  helpers  and 
year  effect  alone  (p=0.0002,  R2=0.09). 

We  have  derived  an  approach  to  predict  recruitment  using  data  that  are  more  easily  collected  at 
broader  spatial  scales  than  age  structure  of  breeders.  The  significance  of  the  year  effect  and 
habitat  attributes  well  outside  of  home  ranges  was  surprising.  More  research  is  required  to 
determine  if  rainfall  and/or  temperature  patterns  might  contribute  to  recruitment.  These  results 
could  lead  to  a  better  SEPM  and  more  effective  habitat  trading  program.  For  example,  some 
parcels  may  not  be  very  valuable  for  RCWs  if  surrounded  by  non-forested  areas  within  3,500  m 
of  cavity  trees.  Ignoring  climatic  variability  could  lead  to  an  under-  or  over-estimate  of  tradable 
credits. 

Overall  Conclusions  and  Recommendations 

The  use  of  Encroachment  Partnering  to  relieve  training  restrictions  due  to  the  presence  of 
endangered  species  highlights  the  importance  of  understanding  how  changing  landscape  patterns 
may  affect  recruitment  and  migration.  Individually-based,  spatially-explicit  population  models 
can  be  used  to  determine  how  changing  habitat  area  and  /  or  connectivity  will  affect  population 
viability.  These  models  allow  interaction  between  rates  of  recruitment  within  patches  with  rates 
of  migration  among  patches.  We  applied  new  techniques  to  improve  how  SEPMs  can  be  used 
for  Encroachment  Partnering.  By  integrating  LEA  and  POM  we  found  that  existing  monitoring 
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data  reduced  enough  uncertainty  regarding  dispersal  to  prioritize  trading  decisions.  However,  we 
also  observed  that  examining  habitat  trades  at  this  spatial  extent  always  resulted  in  increasing 
habitat  fragmentation  effects. 

The  current  study  estimated  habitat  fragmentation  effects  relative  to  the  spatial  distribution  of 
clusters  on  the  Main  Base  that  have  been  planned  to  meet  MCBCL’s  contribution  to  recovery  of 
the  Coastal  North  Carolina  population.  The  Coastal  North  Carolina  population  also  includes 
RCW  clusters  in  Croatan  National  Forest  and  Holly  Shelter  State  Game  Lands,  which  are  located 
within  the  Onslow  Bight  Landscape.  Determining  how  habitat  trading  would  change  levels  of 
habitat  fragmentation  across  the  entire  Onslow  Bight  landscape  is  beyond  the  scope  of  a  one  year 
SEED  grant.  However,  habitat  trades  that  increase  habitat  fragmentation  relative  to  the  current 
population  on  MCBCL  may  actually  decrease  fragmentation  across  the  entire  Coastal  North 
Carolina  recovery  population.  In  order  to  determine  if  DoD  can  spread  RCW  liability  across  the 
Onslow  Bight  landscape  without  jeopardizing  recovery,  we  need  to  understand  how  landscape 
patterns  affect  recruitment  and  migration  at  a  broad  geographic  extent. 

In  order  to  capture  habitat  fragmentation  effects,  LEA  assigns  a  conservation  value  to  trades  at  a 
local  scale  based  on  marginal  changes  in  ecological  services  at  the  landscape  (regional)  scale, 
relative  to  service  levels  expected  in  an  unfragmented  or  baseline  landscape.  We  defined 
baseline  for  MCBCL  as  the  apportionment  of  genetic  variance  within  and  among  breeding 
groups  afforded  by  their  recovery  plan.  By  applying  the  concept  to  the  entire  Onslow  Bight,  we 
might  define  baseline  as  the  apportionment  of  genetic  variance  within  and  among  breeding 
groups,  averaged  across  all  three  management  units  (MUs),  afforded  by  reaching  property- 
specific  recovery  goals  (USFWS  2003a).  In  contrast  to  the  baseline  defined  for  MCBCL,  the 
baseline  for  the  Onslow  Bight  does  not  provide  the  most  clumped  distribution  of  habitat 
possible,  merely  the  most  clumped  distribution  possible  on  public  lands.  We  propose  a 
landscape  genetics  study  in  the  Onslow  Bight  to  detennine  if  the  spatial  distribution  of  RCW 
clusters  specified  by  the  Recovery  Plan  provides  adequate  rates  of  migration  or  if  restoration  of 
private  land  would  increase  population  viability,  justifying  the  use  of  Encroachment  Partnering 
to  provide  regulatory  relief. 

We  propose  a  second  study  to  develop  an  SEPM  for  a  different  threatened  or  endangered  species. 
Prior  to  beginning  research  on  our  SEED  grant,  the  Helmholtz  Centre  for  Environmental 
Research  funded  collaboration  between  Douglas  Bruggeman  and  Thorsten  Wiegand  to  establish 
relationships  between  structural  and  biological  descriptors  of  habitat  fragmentation  (Bruggeman 
et  al.  submitted).  The  research  used  the  RCW  SEPM  and  equilibrium  landscape  conditions  to 
establish  relationships  between  landscape  indices  founded  in  the  fundamental  theory  of  point 
pattern  phenomenon  and  population  genetic  theory.  Results  indicated  that  the  second  order 
spatial  structure  of  habitat  is  well  correlated  to  behavioral-  and  demographic-rules  coded  into  the 
RCW  SEPM.  We  propose  applying  the  approach  “in  reverse”  for  a  new  species.  In  other  words, 
using  data  describing  the  second  moment  approximation  of  alternative  land  cover  types  and 
population  genetic  data  to  begin  constructing  an  SEPM  for  new  species.  We  stress  that 
relationships  between  landscape  indices  and  population  genetic  structure  probably  do  not  provide 
an  adequate  surrogate  for  SEPMs  when  evaluating  management  decisions  that  change  habitat 
spatial  structure  over  time. 
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1.  Phase  I.  Reducing  Uncertainty  in  the  Red-Cockaded  Woodpecker 
Spatially-Explicit  Population  Model  with  Pattern-Oriented  Modeling 

1.1  Objective 

Uncertainty  regarding  critical  biological  processes  required  to  maintain  persistent  subdivided 
populations  should  be  the  focus  of  any  scientifically  defensible  habitat  trading  program.  Habitat 
trading  programs  will  change  how  habitats  are  allocated  for  small  and  sometimes  declining 
populations.  Further,  landscape  patterns  that  may  facilitate  or  restrict  movement  of  individuals 
critical  for  population  persistence  will  be  changing  over  time.  Individually-based  and  spatially- 
explicit  population  models  (SEPMs)  are  one  approach  to  test  the  response  of  populations  to  non¬ 
equilibrium  conditions,  created  by  habitat  trading  programs.  The  objective  of  this  phase  was  to 
develop  an  uncertainty  analysis  for  the  Red-Cockaded  Woodpecker  (RCW)  SEPM  currently 
used  at  Marine  Corps  Base  Camp  Lejeune  (MCBCL). 

1.2  Background 

A  SEPM  is  comprised  of  linked  submodels  that  simulate  how  different  biological  processes  are 
affected  by  land  cover  patterns.  Initially,  SEPMs  were  heavily  criticized  for  their  complexity, 
large  data  requirements,  and  reliance  on  uncertain  parameters  (Beissinger  and  Westphal  1998). 
However,  Pattern-Oriented  Modeling  (POM)  has  been  derived  to  construct  and  validate  these 
models  with  less  data  than  previously  imagined  (Grimm  et  al.  2005;  Wiegand  et  al.  2004).  POM 
accomplishes  this  by  using  additional  information  encoded  in  observed  population-level  patterns 
for  indirect  parameter  estimation.  A  pattern  is  defined  as  anything  beyond  random  variation,  or 
any  signal  beyond  noise,  and  is  used  as  an  indicator  of  the  system’s  internal  behavior  (Grimm  et 
al.  1996).  Pattern-Oriented  Modeling  also  is  a  method  of  comparing  the  ability  of  different 
SEPMs  to  reproduce  biological  patterns  observed  in  nature  (Grimm  et  al.  2005). 

While  POM  can  be  used  to  derive  the  structure  of  SEPMs,  or  the  types  of  submodels  included 
and  how  they  interrelate  (Grimm  et  al.  1996),  we  focused  on  using  POM  to  reduce  parameter 
error  (Wiegand  et  al.  2004).  Parameter  error  is  the  uncertainty  associated  with  the  parameter 
values  used  in  a  submodel.  POM  can  be  used  to  reduce  parameter  error  by  comparing  the  ability 
of  alternative  parameterizations  of  a  model  to  reproduce  patterns  observed  in  nature.  First,  one 
must  detennine  the  range  of  possible  parameter  values  based  on  expert  knowledge  or  parameters 
derived  in  landscapes  where  data  are  more  abundant.  Second,  the  range  of  possible  parameter 
values  is  subdivided  into  equidistant  values.  Third,  stratified  sampling  of  possible  parameter  sets 
without  replacement  is  used  to  combine  all  possible  parameter  values  that  are  treated  as 
uncertain.  Fourth,  the  SEPM  is  run  with  several  thousand  different  parameterizations  to  ensure 
exhaustive  sampling  of  possible  parameter  combinations.  Fifth,  parameterizations  unable  to 
reproduce  RCW  patterns  observed  on  MCBCL  are  removed  from  consideration  through 
statistical  comparisons  between  patterns  produced  by  the  simulation  (i.e.,  expected  patterns)  and 
patterns  observed  on  MCBCL  (i.e.,  observed  patterns).  This  process  will  be  referred  to  as 
“filtering”  parameterizations. 

We  applied  POM  to  the  RCW  SEPM  originally  described  by  Letcher  et  al.  (1998).  The  Letcher 
model  assumes  the  population  consists  of  breeders,  helpers,  floaters,  and  fledglings  (Fig.  1.1). 
The  probability  of  an  individual  transitioning  among  these  classes  depends  on  demographic, 
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behavioral,  and  landscape  factors.  Our  study  focuses  on  uncertainty  regarding  how  dispersal 
behaviors  affect  population  processes.  The  Letcher  SEPM  assumed  that  when  birds  leave  their 
natal  territory  they  choose  a  direction  of  travel  at  random  and  then  travel  in  a  straight  line. 
However,  field  biologists  have  observed  that  many  fledglings  embark  on  forays,  seemingly 
investigating  nearby  breeding  habitat  (Jeffrey  Walters,  personal  communication).  Therefore,  it 
may  be  best  to  assume  that  birds  use  knowledge  of  forest  structure  when  deciding  which 
direction  to  travel.  Others  have  observed  that  birds  often  avoid  crossing  non-forested  areas 
during  foraging  (Connor  and  Rudolph,  1991;  Belisle  and  Desrochers  2002),  seemingly  to 
increase  the  chance  of  finding  breeding  and  foraging  habitat  or  refuge  from  predation. 
Therefore,  birds  may  chose  not  to  travel  straight  in  order  to  avoid  crossing  non-forested  areas. 


Figure  1,1.  Breeding  group  structure  and  demographic  transitions  for  the  Red-Cockaded 
Woodpecker.  Black  boxes  denote  breeding  territories  and  the  grey  box  denotes  the  matrix 
crossed  during  floating  behaviors. 

Ecologists  have  derived  many  ways  of  simulating  dispersal  behaviors.  The  most  simple  is  a 
random  walk  in  which  individuals  choose  each  step  at  random.  Under  a  correlated  random  walk, 
the  turning  angles  between  successive  moves  are  not  statistically  independent  (Zollner  and  Lima 
1999).  A  random  component  remains  in  this  model  because  the  direction  of  travel  selected  is 
based  on  comparing  probabilities  of  turning  (do°,  d45°,  etc...)  to  a  random  number,  often  creating 
a  zigzag  pattern  of  dispersal  even  if  the  probability  of  not  turning  (d0°)  is  large.  To  incorporate 
the  influence  of  forest  cover,  we  will  derive  a  Habitat-based  Correlated  Random  Walk  (HCRW) 
similar  to  that  used  by  Kramer-Schadt  et  al.  (2004).  HCRW  combines  the  probabilities  of 
turning  different  angles  with  the  probability  that  a  bird  will  select  different  land  cover  types  for 
dispersal  (Ricketts  et  al.  2001).  Some  land  cover  types  may  limit  or  facilitate  dispersal  among 
habitat  patches,  thus  varying  in  levels  of  “landscape  friction”,  modeled  as  parameter  “m”.  By 
introducing  a  third  parameter  “B”  that  weighs  the  level  of  importance  birds  place  in  directional 
or  habitat  rules,  the  HCRW  can  be  represented  as: 

P[Cit+ 1  | cJt ] - (1  - B)d'  +  B^m'  / 2> j  .  (LI) 

Where  PfC't+ilC^]  is  the  probability  the  bird  will  be  in  cell  i  at  time  t+1  given  that  it  is  currently 
in  cell  j  now,  d1  is  the  probability  of  turning  toward  cell  i  relative  to  the  direction  of  travel  in  the 
previous  time  step,  and  m1  is  the  landscape  friction  assigned  to  cell  i,  which  is  normalized  by  the 
sum  of  the  friction  values  for  all  8  cells  in  the  floaters  immediate  neighborhood. 
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We  are  uncertain  what  the  appropriate  parameter  values  for  d,  m,  and  B  in  the  HCRW  should  be. 
The  objective  of  Phase  1  was  to  use  POM  to  estimate  parameter  values  for  the  HCRW  using 
RCW  patterns  observed  on  MCBCL.  The  Letcher  model  was  derived  from  monitoring  data  from 
Fort  Bragg  not  Camp  Lejeune.  So  we  examined  uncertainty  in  other  parameters  associated  with 
dispersal,  including  the  number  of  steps  taken,  survival  of  floaters,  and  distance  at  which  a  bird 
can  detect  a  breeding  vacancy.  A  total  of  2 1  parameters  are  treated  as  uncertain.  The  uncertain 
parameter  values  were  filtered  through  multiple  RCW  patterns  observed  on  MCBCL.  Examples 
of  patterns  include  temporal  patterns  of  population  dynamics  and  spatial  patterns  of  territory 
connectivity  and  pairwise  genetic  distances  among  groups.  To  the  best  of  our  knowledge  this 
study  represents  the  first  application  of  population  genetic  patterns  to  POM.  Below  we  describe 
POM  in  the  context  of  calibrating  the  RCW  SEPM  to  the  patterns  observed  on  Camp  Lejeune. 

1.3  Methods 

1.3.1  Landscape  Map 

The  Camp  Lejeune  Land  Cover  Map  was  used  for  the  simulation.  The  Land  Cover  Map 
provided  a  general  description  of  forest  cover  updated  every  ten  years,  and  information  on  the 
characteristics  of  non-forested  areas  thought  to  create  landscape  friction.  The  Land  Cover  Map 
classifies  areas  on  MCBCL  into  13  types;  we  condensed  this  list  down  to  7  types  by  combining 
similar  areas  (Table  1.1).  Land  Cover  maps  for  offsite  areas  were  not  available  so  we  assumed 
that  these  areas  are  homogeneous  with  regard  to  RCW  dispersal  behaviors  and  designated  them 
as  an  eighth  type. 

Table  1.1.  Modified  classification  of  Camp  Lejeune  Land  Cover  Map  used  to  estimate  friction  of 
land  cover  types  for  dispersal. _ 


Code 

Camp  Lejeune  Land  Cover  Classes 

1 

Predominantly  Pine  Forest 

2 

Mixed  Pine  and  Hardwood  Forest 

3 

Bottomland  Hardwood  /  Upland  Hardwood  Forest 

4 

Forest  plantation  under  10  years 

5 

Bare  ground  /  cropland  /  grasses  /  marshland  /  scrub,  shrub 

6 

River,  Stream,  Lake 

7 

Business  or  Commercial  area  /  Residential  Area 

8 

Off  site  areas 

Combined  classes  designated  by 


The  Land  Cover  shape  file  was  converted  into  a  raster  file  in  ArcView  3.3  using  Spatial  Analyst, 
for  use  in  dynamic  simulation  modeling.  A  100  x  100  m  cell  size  (1  ha)  was  chosen  as  the  grain 
for  the  simulation  study.  This  value  represents  a  tradeoff  between  a  need  to  limit  computer 
simulation  time  to  test  many  uncertain  parameters  (i.e.,  summarizing  more  area  with  fewer  cells) 
and  level  geographic  complexity  perceived  by  RCWs  when  making  dispersal  choices.  The 
Letcher  model  assumes  that  birds  can  detect  breeding  vacancies  within  3  km,  so  we  felt  a  1  ha 
cell  size  was  sufficiently  conservative.  MCBCL  also  provided  a  shape  file  with  approximate 
cluster  centers  for  all  existing  RCW  territories.  After  creating  the  land  cover  raster  file,  we  used 
ArcView  3.3  Spatial  Analyst  to  generate  a  raster  file  of  cluster  centers,  each  consisting  of  a  1-ha 
cell. 
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1 .  Pine  Forest 


2.  Pine  /  Hardwood  Forest 

3.  Hardwood  Forest 

4.  Forest  <10  years  old 


5.  Non-forested 

6.  Water 


7.  Developed  Areas 

8.  Off-site  Areas 


Figure  1.2.  Camp  Lejeune  land  cover  map  used  in  the  spatially-explicit  population  model.  Cell  size  is  100  m  x  100  in,  or  1- 
hectare.  White  dots  designate  approximate  cluster  centers.  All  clusters  plotted  were  available  for  RCWs  as  of  2005. 
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1.3.2  Initial  Conditions 

Landscape  pattern  data  are  updated  at  least  every  10  years  in  MCBCL’s  GIS  system.  Therefore, 
we  calibrated  the  RCW  SEPM  to  landscape  patterns  for  the  years  1995-2005.  We  used  the  RCW 
monitoring  data,  provided  thanks  to  Dr.  Jeffrey  Walters,  VPI&SU,  to  determine  the  population 
structure  present  at  year  1995.  We  identified  each  bird’s  natal  territory,  past  and  current 
breeding  territory,  present  location,  age,  status  (fledgling,  helper,  floater  or  breeder),  and  sex. 
There  were  192  birds  in  the  initial  population.  All  immigrants  to  Lejeune  between  1995  and 
2005  were  also  identified  -  65  in  total.  The  location,  age,  and  sex  were  unknown  for  some 
floaters  and  immigrants,  so  we  assigned  these  at  random  before  the  start  of  each  simulation, 
which  was  replicated  50  times.  We  felt  this  is  an  effective  way  to  remove  any  bias  caused  by 
randomly  assigning  unknown  characteristics  to  birds  because  we  summarize  the  results  of  one 
SEPM  with  an  average  of  50  replicate  simulations. 

In  order  to  estimate  the  population  genetic  structure  of  RCWs  on  MCBCL  we  used  an  infinite 
alleles  model  (Lacy  and  Lindenmayer,  1995).  The  infinite  alleles  model  assumes  that  every 
founding  individual  is  heterozygous  and  contains  two  unique  alleles  (i.e.,  total  alleles  =  2  x 
number  of  breeders  in  the  founding  population).  By  making  this  assumption  we  can  examine 
how  the  processes  of  inbreeding,  genetic  drift  and  gene  flow  interact  in  a  spatially  subdivided 
population.  To  apply  the  infinite  alleles  model,  we  identified  all  founding  individuals  or 
individuals  observed  between  1986  and  1995  that  reproduced  but  for  which  their  parentage  is 
unknown.  Each  was  assigned  2  unique  alleles  at  a  locus,  therefore,  founders  were  assumed  to  be 
unrelated.  Next,  we  applied  a  gene  drop  which  simulates  the  transmission  of  alleles  from  parents 
to  offspring  assuming  Mendelian  inheritance  (i.e.,  their  offspring  had  an  equal  probability  of 
inheriting  each  of  the  two  alleles).  This  stochastic  process  was  repeated  for  the  years  1986  to 
1 995  at  the  start  of  each  simulation. 

The  internal  structure  of  the  model  must  closely  approximate  conditions  and  processes  observed 
in  nature  to  apply  POM  (Wiegand  et  al.  2004).  Therefore,  all  immigrants  to  MCBCL  between 
1995  and  2005  were  also  assigned  2  unique  alleles  at  a  locus.  In  other  words,  we  treated  all 
observed  immigrants  as  “founders”  as  it  is  possible  that  they  may  successfully  reproduce  some 
time  during  the  simulation.  Further,  when  recruitment  clusters  were  made  available  under  the 
MCBCL  RCW  Management  Plan  or  when  budding  or  pioneering  behaviors  led  to  new  groups, 
the  territories  were  included  in  the  simulation  at  the  appropriate  time. 

1.3.3  Parameters  Used  in  the  Simulation 

1.3. 3.1  Recruitment  and  Survival.  The  Letcher  model  uses  an  age-based 
recruitment  submodel  in  which  the  probability  that  birds  attempt  to  nest,  probability  of 
nest  success,  and  number  of  fledglings  produced  is  a  function  of  the  age  of  the  male  and 
female  breeder  and  the  number  of  helpers  present.  We  did  not  have  data  regarding  which 
groups  attempted  to  nest  and  of  those  how  many  were  successful.  Therefore,  we  used  the 
point  estimates  for  probability  that  birds  attempt  to  nest  and  probability  of  nest  success 
from  the  Sandhills  population  used  previously  in  the  Letcher  model.  However,  we  did 
directly  estimate  the  number  of  fledglings  produced  from  the  Camp  Lejeune  population 
using  the  regression  equation  #3  from  Letcher  et  al.  (1998),  for  the  years  1995-2005: 

FL  =  64  +  b$e  ^Sem  +  b^e  +b-jHP .  (eqn  1.2) 
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The  resultant  parameter  estimates  are  reported  in  Table  1.2  (F(3,458)=14.65,  p<0. 00001). 

Estimates  of  the  probability  of  individual  survival  on  Camp  Lejeune  were  recently 
estimated  by  Walters  (2004)  for  the  years  1986-2000  (Table  1.2).  We  used  these 
estimates  for  our  simulation  study.  However,  estimates  for  female  helpers  and  floaters  of 
both  sexes  were  not  included  in  the  previous  study.  The  presence  of  female  helpers  is  a 
recent  observation  and  may  be  due  to  population  growth  rates  exceeding  rates  at  which 
recruitment  clusters  have  been  provisioned  on  MCBCL  (Walters  2004).  This  change  in 
population  structure  was  not  the  focus  of  our  research  so  we  assumed  that  female  helpers 
have  the  same  survival  probability  as  male  helpers,  because  floaters  of  both  sexes  have 
also  been  shown  to  have  equal  survival  probability  (Letcher  et  al.  1998).  Probability  of 
floater  survival  for  both  sexes  was  treated  as  uncertain  because  we  are  concerned  about 
understanding  floater  behaviors  in  this  landscape.  We  included  a  range  of  values  for 
probability  of  seasonal  survival,  0.8  to  1.0,  because  the  average  observed  in  the 
Sandhills’  population  was  0.89  (Letcher  et  al.  1998). 

1.3. 3.2  Dispersal  Behaviors.  The  RCW  SEPM  simulates  dispersal  behaviors  in 
two  steps.  The  first  step  is  to  detennine  if  individuals  delay  dispersal  and  stay  as  helpers. 
To  simulate  this  step  we  used  point  estimates,  as  was  done  for  the  Letcher  model,  from 
the  years  1995-2005  (Table  1.3)  because  we  were  most  interested  in  how  land  use 
mediates  dispersal.  We  did  make  an  important  structural  improvement  by  allowing 
females  to  delay  dispersal.  Ignoring  this  structural  difference  between  the  model  and 
conditions  observed  on  Camp  Lejeune  would  prevent  POM  from  effectively  reducing 
uncertainty. 

The  second  step  is  simulating  floating  or  dispersal  behaviors  in  a  heterogeneous 
landscape.  We  must  first  determine  the  number  of  steps  taken  each  season.  The  Letcher 
model  assumes  that  females  and  first  year  male  floaters  can  disperse  farther  than  older 
male  floaters.  We  treated  these  values  as  uncertain,  allowing  dispersal  steps  to  vary  from 
1-10  km  per  season  for  females  and  first  year  male  floaters,  and  1-8  km  per  season  for  the 
older  male  floaters  (Table  1.3). 

Next,  we  need  to  simulate  where  each  bird  moves  to  with  each  step,  using  the  HCRW 
model  described  above.  We  are  uncertain  regarding  the  preferences  RCW  have  for 
different  land  cover  types  during  dispersal,  or  conversely  the  level  of  friction  land  cover 
types  represent.  Therefore,  we  used  a  range  of  values  for  the  degree  of  landscape  friction 
the  eight  land  cover  types  (Table  1.1)  represent.  The  values  for  “m”  (eqn  1.1)  range  from 
0.01,  meaning  that  birds  have  a  low  probability  of  entering,  to  1.0.  Similarly,  we  are 
uncertain  regarding  the  turning  decisions  birds  make  during  dispersal.  From  a  human’s 
perspective  on  the  ground,  RCWs  are  often  seen  to  be  dispersing  in  a  straight  line 
(Letcher  et  al.,  1998).  However,  a  human’s  perceptual  distance,  without  the  aide  of  radio 
telemetry,  is  limited,  so  directional  preferences  should  be  treated  as  uncertain.  When 
assigning  values  to  all  possible  “d”  parameters,  we  must  recognize  that  these  values  are 
not  statistically  independent,  when  turning  a  bird  may  turn  to  the  left  or  the  right,  and 
they  should  sum  to  unity  (do°  +  2(d45°)  +  2(d9o°)  +  2(d 1 35°)  =  1).  Therefore,  we  started  by 
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randomly  assigning  the  value  for  do°,  tendency  for  straight  movement.  Then  a  random 
number  between  do  and  1  was  selected  and  divided  by  2  to  determine  d45°.  This  step  was 
then  repeated  for  the  remaining  two  possible  directions. 

1.3. 3. 3  Competition  for  Breeding  Vacancies.  Competition  for  breeding 
vacancies  is  an  important  biological  process  that  ultimately  detennines  if  a  dispersing 
bird  will  contribute  to  population  growth.  The  competition  rules  used  by  the  Letcher 
model  were  reproduced  here.  For  males,  it  assumes  that  the  oldest  male  helpers  will 
always  win  the  breeding  vacancy  in  the  territory  in  which  they  are  helping.  If  no  helpers 
are  present,  then  helpers  and  floaters  within  a  specified  search  radius  will  compete  for  the 
vacancy.  The  model  assumes  that  the  closest  individual  to  the  vacant  territory  will  win, 
but  if  individuals  are  equi-distant,  the  oldest  individual  wins. 

For  females,  the  Letcher  model  assumes  that  the  oldest  individual  within  a  specified 
search  radius  will  win.  However,  by  adding  female  helpers  we  had  to  examine  how 
competition  rules  could  be  extended  to  a  new  class  of  birds.  Investigating  the  RCW 
monitoring  data  revealed  that  unlike  male  helpers,  female  helpers  rarely  inherit  their  natal 
territory.  This  occurred  twice,  on  Cluster  14  in  1990  and  on  Cluster  28  in  1989.  Both 
females  were  born  in  1987,  only  the  second  year  of  bird  banding  studies  on  Camp 
Lejeune,  and  their  parentage  was  uncertain.  On  cluster  14  the  female  helper  inherited  the 
position  despite  there  being  a  male  fledgling  present,  which  stayed  as  a  helper  in  that 
territory.  In  territory  28,  only  female  fledglings  were  produced  prior  to  the  occurrence  of 
a  breeding  vacancy.  In  both  cases,  a  new  male  breeder  was  observed.  In  2003  at  Cluster 
48,  a  female  fledgling  assumed  the  breeding  position  after  her  parents  and  the  helper 
disappeared  from  the  cluster.  These  events  were  rare  and  may  reflect  predation  events 
that  are  not  well  characterized  by  the  model.  Therefore,  we  assumed  that  female 
fledglings,  helpers,  and  floaters  could  not  win  a  vacancy  in  their  natal  cluster.  Further,  as 
assumed  in  Letcher  et  al.  (1998),  we  assumed  that  when  a  male  helper  takes  a  breeding 
position,  the  female  breeder  dispersers  and  can  not  re -win  a  vacancy  where  she 
successfully  bred  in  the  past. 

The  ability  of  birds  to  detect  and  acquire  breeding  vacancies  will  have  a  large  impact  on 
the  persistence  and  population  structure  in  a  fragmented  landscape.  In  absence  of 
empirical  estimates  of  a  bird’s  search  radius  the  Letcher  model  assumes  that  all 
fledglings,  helpers,  and  floaters  of  both  sexes  can  compete  for  breeding  vacancies  within 
3  km  of  their  current  location.  Such  a  large  search  radius  may  be  plausible  for  female 
floaters  and  helpers  because  it  is  believed  that  solitary  males  drum  on  pine  trees  to  attract 
distant  females  (Conner  et  al.,  2001).  Also,  male  and  female  helpers  are  known  to 
conduct  forays.  However,  male  floaters  investigating  new  areas  may  have  a  very  small 
search  radius.  We  used  POM  to  test  the  assumption  that  floaters  and  helpers  of  both 
sexes  have  a  similar  search  radius.  Therefore,  we  expanded  the  model  to  include  4  search 
radii  (one  for  each  sex  for  helpers  and  floaters)  and  varied  the  values  from  0.1  to  6  km 
(Table  1.3).  We  assumed  that  birds  are  able  to  detect  vacancies  at  this  distance  are 
willing  to  cross  all  land  cover  types  to  acquire  a  breeding  vacancy  (e.g.,  Norris  and 
Stutchbury  2001). 
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Table  1.2.  Recruitment  and  survival  parameter  values  used  in  Red-Cockaded  Woodpecker 
spatially-explicit  population  model. _ _ _ 


Parameter 

Value 

Source 

a  Female  nesting  attempt 

1.14569 

Letcher  et  al.  1998 

bo  Nest  success  intercept 

2.3404 

Letcher  et  al.  1998 

bi  Nest  success  male’s  effect 

-2.4295 

Letcher  et  al.  1998 

bi  Nest  success  female’s  effect 

-1.1527 

Letcher  et  al.  1998 

b3  Nest  success  helper’s  effect 

-1.1012 

Letcher  et  al.  1998 

b4  Mean  fledgling  intercept 

1.8661 

This  study 

bs  Mean  fledgling  male’s  effect 

-1.4351 

This  study 

b6  Mean  fledgling  female’s  effect 

-0.3288 

This  study 

by  Mean  fledgling  helper’s  effect 

0.1766 

This  study 

Seasonal  probability  of  survival  -  Male  Fledgling 

0.88 

Walters  2004 

Seasonal  probability  of  survival  -  Male  Helper 

0.94 

Walters  2004 

Seasonal  probability  of  survival  -  Male  Breeder 

0.95 

Walters  2004 

Seasonal  probability  of  survival  -  Male  Floater 

p 

i 

00 

o 

This  study 

Seasonal  probability  of  survival  -  Female  Fledgling 

0.81 

Walters  2004 

Seasonal  probability  of  survival  -  Female  Helper 

0.94 

Walters  2004 

Seasonal  probability  of  survival  -  Female  Breeder 

0.94 

Walters  2004 

Seasonal  probability  of  survival  -  Female  Floater 

o 

l 

00 

o 

This  study 

Maximum  age 

17 

Conner  et  al.  2001 
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Table  1.3.  Dispersal  and  competition  parameter  values  used  in  Red-Cockaded  Woodpecker 
spatially-explicit  population  model. _ _ _ 


Parameter 

Value 

Source 

Probability  of  Male  Natal  Dispersal 

0.568 

This  study 

Probability  of  Female  Natal  Dispersal 

0.874 

This  study 

Female  floater  and  Male  natal  dispersal  speed  (per  season) 

1  -  10  km 

This  study 

Male  floater  dispersal  speed  (per  season) 

1  -  8  km 

This  study 

ml  Probability  of  entering  predominantly  pine  forest 

0.01  -  1.0 

This  study 

m2  Probability  of  entering  mixed  pine  and  hardwood  forest 

0.01  -  1.0 

This  study 

m3  Probability  of  entering  bottomland  hardwood  /  upland 
hardwood  forest 

0.01  -  1.0 

This  study 

m4  Probability  of  entering  forest  plantation  under  10  years 

0.01  -  1.0 

This  study 

m5  Probability  of  entering  bare  ground  /  cropland  /  grasses  / 
marshland  /  scrub,  shrub 

0.01  -  1.0 

This  study 

m6  Probability  of  entering  river,  stream,  lake 

0.01  -  1.0 

This  study 

m7  Probability  of  entering  business  or  commercial  area  / 
residential  area 

0.01  -  1.0 

This  study 

m8  Probability  of  entering  off  site  areas 

0.01  -  1.0 

This  study 

do  Probability  of  traveling  in  same  direction  that  was  chosen 
in  last  time  step 

0.01  -  1.0 

This  study 

d45°  Probability  of  turning  45°  away  from  direction  traveled  in 
previous  time  step 

0.01  -  1.0 

This  study 

dgo°  Probability  of  turning  90°  away  from  direction  traveled  in 
previous  time  step 

0.01  -  1.0 

This  study 

di35°  Probability  of  turning  135°  away  from  direction  traveled 
in  previous  time  step 

0.01  -  1.0 

This  study 

B  Strength  of  preference  for  habitat-based  or  directional 
dispersal  rules. 

B  =  1 ,  only  land  cover  type  matters 

B  =  0,  only  directional  rules  matter 

0-1 

This  study 

Perceptual  distance  of  Male  Helpers 

0.1-6  km 

This  study 

Perceptual  distance  of  Male  Floaters 

0.1-6  km 

This  study 

Perceptual  distance  of  Female  Helpers 

0.1-6  km 

This  study 

Perceptual  distance  of  Female  Floaters 

0.1-6  km 

This  study 
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1.3.4  Simulation  Conditions 

The  SEPM  was  run  with  100,000  different  parameterizations,  generated  by  a  Latin  hypercube,  to 
ensure  exhaustive  sampling  of  possible  parameter  combinations.  A  Latin  hypercube  design  is 
used  to  provide  a  stratified  sampling  of  possible  parameter  sets  without  replacement  (McKay  et 
al.  2000).  The  range  of  values  for  each  of  the  21  uncertain  parameters  reported  above  were 
subdivided  into  2 1  equidistant  values.  The  Latin  hypercube  was  then  used  to  randomly  assemble 
parameter  values  for  the  21  uncertain  parameters.  This  represented  one  possible 
parameterization.  Each  parameterization  was  then  simulated  for  1 1  years  (1995-2005)  and 
replicated  50  times. 


1.3.5  Filtering:  Model  Calibration  by  Comparing  Observed  and  Expected  Patterns 

The  RCW  SEPM  was  calibrated  to  MCBCL  by  filtering  the  100,000  parameterizations  based  on 
their  ability  to  reproduce  patterns  observed  in  nature.  Pattern  fulfillment  indices  are  developed 
that  measure  the  deviance  between  observed  and  simulated  patterns.  Multiple  indices  were  used 
because  compensatory  mechanisms  such  as  low  probability  of  floater  survival  being  offset  by 
high  dispersal  speed  may  result  in  several  different  parameterizations  being  able  to  reproduce  the 
same  pattern.  Therefore,  observed  patterns  were  chosen  to  reflect  multiple  biological  processes 
to  remove  the  effects  of  compensatory  mechanisms.  Further,  we  must  derive  pattern  fulfillment 
criteria  that  recognize  that  the  patterns  observed  may  not  reflect  underlying  biological  process 
but  instead  may  be  due  to  random  chance.  Pattern  fulfillment  criteria  are  the  statistical  criteria 
used  to  remove  parameterizations  that  simulate  patterns  that  do  not  differ  significantly  from 
patterns  resulting  from  random  processes.  The  pattern  fulfillment  criteria  represent  the 
“porosity”  used  during  filtering.  Filtering  is  then  the  removal  of  parameterizations  with  pattern 
fulfillment  indices  not  able  to  satisfy  the  pattern  fulfillment  criteria 


1.3. 5.1  Index  1.  Potential  Breeding  Groups  (PBG)  serve  as  the  regulatory 
benchmark  for  RCW  Recovery  (USFWS  2003a).  We  tested  the  ability  of  different 
parameterizations  to  reproduce  the  time  series  of  PBG  observed  on  MCMCL  from  1995 
to  2005  (Fig.  1.3).  We  measured  the  deviance  between  the  PBGs  present  at  the  start  of 
the  breeding  season  in  the  simulation  (S[PBGt])  and  those  observed  on  MCBCL 
(0[PBGt])  for  each  of  the  1 1  years.  An  error  measure,  EPBG,  was  derived  using  the  root 
squared  mean  deviation  between  the  observed  and  expected  time  series  (Wiegand  et  al. 
2004): 


ePBG=.±-  2 1 tcAPBG, ] - S[PBGt ]) 

V 11 (=1995 


(eqn.  1.3) 


Because  each  parameterization  was  replicated  50  times,  the  error  value  for  each  replicate 
was  estimated.  The  average  EPBG  from  these  50  replicates  was  used  to  characterize  the 
suitability  of  the  parameterization  (Wiegand  et  al.  2004). 


The  pattern  fulfillment  criteria,  or  filter  porosity,  for  EPBG  was  determined  using  a 
randomization  approach.  Random  time  series  of  PBGs  were  assembled  by  randomly 
assigning  a  year  between  1995  and  2005  to  each  of  the  1 1  observed  PBGs  (i.e.,  random 
pennutation),  generating  random  time  series,  R[PBGt].  EPBGRand  was  then  calculated: 
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E PBG  = 
Rand 


™{c{PBGt]-R[PBGt]) 

/— 1995 


(eqn.  1.4) 


To  ensure  consistency  in  analysis  between  the  simulated  and  random  time  series,  random 
pennutation  of  the  time  series  was  repeated  50  times  and  the  final  EPBGRand  represents  an 
average  of  the  50.  This  procedure  was  repeated  10,000  times  to  produce  a  distribution  of 
EPBGRand-  The  pattern  fulfillment  criteria  for  PBG  specified  that  EPBG  must  be  less  than 
95%  of  all  EPBGRand  or  a  root  squared  mean  deviation  of  15.2.  This  very  strict  criterion 
was  chosen  because  PBGs  have  regulatory  significance  and  sampling  of  potential 
breeding  habitat  for  RCWs  is  conducted  thoroughly  each  year.  Therefore,  we  should 
expect  our  SEPM  to  predict  this  trend  with  little  error. 


1.3. 5.2  Index  2.  The  number  of  male  and  female  floaters  observed  each  year  on  MCBCL 
served  as  our  second  index  (Fig.  1.3).  Index  2  was  developed  exactly  as  Index  1  because 
they  are  both  time  series  data.  Criteria  for  males,  for  example,  were  derived  by 
substituting  0[mFloatert],  S[mFloatert],  and  R[mFloatert]  in  the  equations  to  derive 

EmFloater j  T-'mFloater 

and  E  Rand. 


Unlike  PBG,  we  expect  a  large  amount  of  observational  error  in  OfFloater]  on  MCBCL, 
because  the  monitoring  program  focuses  on  individuals  located  in  clusters,  not  floaters. 
However,  it  is  important  to  include  this  pattern  when  filtering,  because  parameterizations 
including  high  floater  survival  may  lead  to  unrealistic  rates  of  immigration  and  territory 
occupancy.  To  account  for  the  large  amount  of  observational  error,  the  pattern 
fulfillment  criteria  for  Floaters  specified  that  El  loalei  be  within  the  range  of  error  values 
generated  when  the  time  series  was  randomized.  In  other  words,  EFloatei  must  be  no 
greater  than  the  maximum  error  observed  when  the  time  series  were  randomized, 
max(EmFloaterRand)  =  2.9  and  max(EffloaterRand)  =  7.95. 
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Figure  1.3.  Observed  temporal  patterns  on  MCBCL  used  to  establish  pattern 
fulfillment  indices  for  Potential  Breeding  Groups  and  Floaters. 
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1.3. 5.3  Index  3.  Group  size  observed  at  each  territory  and  each  year  was  also 
used  as  a  pattern.  We  defined  group  size  as  the  number  of  breeding  adults  and  helpers  in 
a  territory.  To  test  the  ability  of  one  parameterization  to  reproduce  the  temporal  and 
spatial  distribution  of  breeding  adults  and  helpers  observed  on  MCBCL,  a  Pearson 
correlation  coefficient  was  calculated  by  comparing  group  size  observed  at  each  territory 
and  each  year  in  the  simulation  with  those  observed  in  the  field.  We  will  refer  to  this  as 
our  “observed  test  statistic”,  robs-  An  average  correlation  coefficient  from  the  50 
replicates  was  used  to  characterize  the  ability  of  each  parameterization  to  reproduce 
patterns  observed  on  MCBCL. 

To  derive  the  pattern  fulfillment  criterion,  we  calculated  a  distribution  of  test  statistics  by 
randomly  exchanging  group  size  among  years  and  territories  (i.e.,  random  permutation) 
and  correlating  this  random  pattern  with  the  observed  pattern,  rRand.  Recall  that  the 
simulation  introduces  new  territories  to  the  landscape  to  reflect  when  recruitment  clusters 
were  installed  by  MCBCL  and  when  new  clusters  started  by  RCWs,  through  budding  and 
pioneering  behaviors,  were  observed.  Therefore,  to  control  for  this  autocorrelation 
between  the  simulation  and  field  conditions,  values  for  group  size  were  randomly 
exchanged  only  among  territories  that  were  available,  whether  occupied  or  not.  The 
random  permutation  was  repeated  10,000  times  and  rRand  was  calculated  for  each 
pennutation.  The  null  hypothesis  tested  was  that  the  correlation  between  simulated  and 
observed  distribution  of  group  sizes  is  no  different  than  would  be  due  to  random 
processes.  We  can  calculate  the  probability  that  robs  resulted  due  to  random  processes 
using  the  following  equation  (Manly  1986): 


NrRand  ^  rObs 

10,000 


(eqn  1.5) 


Where,  NrRand>  r0bs  is  the  number  of  times  rRand  was  observed  to  be  greater  than  or  equal 
to  r0bs-  The  maximum  value  for  rRand  was  0.56.  Therefore,  the  probability  of  observing 
robs  equal  to  0.56  resulting  from  random  spatial  processes  is  1/10,000  or  p=0.0001.  The 
pattern  fulfillment  criterion  for  group  size  specified  that  parameterizations  must  produce 
robs  greater  than  0.56  (p<0.0001). 

1.3. 5.4  Index  4.  Connectivity  of  territories  was  used  as  the  fourth  pattern. 
Connectivity  was  estimated  only  for  fledglings,  helpers,  and  floaters  that  transitioned  into 
a  breeding  position  in  both  the  simulation  and  observed  demographic  transitions  on 
MCBCL.  Connectivity  was  estimated  for  males  and  females  separately,  simply  as  the 
number  of  individuals  exchanged  among  territories  between  1995  and  2005  (e.g.,  xy, 
where  x  =  number  of  females  originating  at  territory  i  and  achieving  breeding  status  in 
territory  j).  MCBCL  has  numbered  each  RCW  cluster  from  1  to  90.  Therefore,  a  square 
90x90  matrix  of  connectivity  values  was  used  to  summarize  pair  wise,  directional 
connectivity  among  territories.  A  diagonal  value  in  the  matrix  for  males  then  indicates 
the  number  of  times  a  male  helper  inherited  the  territory.  A  Mantel  test,  which  estimates 
correlation  among  matrices,  was  used  to  test  for  the  degree  of  similarity  between  the 
observed  connectivity  matrix  and  expected  connectivity  matrix  generated  in  the 
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simulation  (Bailey  and  Gatrell  1995).  An  average  Mantel  correlation  coefficient  from  the 
50  replicates  was  used  to  characterize  the  ability  of  each  parameterization  to  reproduce 
patterns  observed  on  MCBCL,  MrCobs- 

To  derive  the  pattern  fulfillment  criterion,  we  calculated  a  distribution  of  Mantel  test 

statistics  by  randomly  exchanging  connectivity  values  among  territories  (i.e.,  random 

c 

pennutation)  and  correlating  this  random  pattern  with  the  observed  pattern,  Mr  Rand.  This 
was  performed  10,000  times  as  explained  above  for  Group  Size.  The  pattern  fulfillment 
criterion  for  connectivity  specified  that  parameterizations  must  produce  Mr  obs  greater 
than  0.1 16  (p<0.0001)  for  males  and  0.065  (p<0.0001)  for  females. 

1.3. 5. 5  Index  5.  Pair  wise  minimum  genetic  distance  (D,j)  among  all  territories 
was  also  used  as  a  pattern  (Nei  1973).  A  mathematical  definition  is  given  in  equation  A4 
in  Appendix  A.  When  D,,  is  close  to  zero,  then  alleles  are  present  in  territory  i  and  j  at 
very  similar  frequencies,  indicating  that  gene  flow  between  the  groups  is  high.  Recall 
that  an  infinite  alleles  model  of  genetic  variance  was  included  in  each  simulation,  and 
that  a  gene  drop  was  used  to  simulate  the  transmission  of  alleles  from  parents  to 
offspring.  Therefore,  a  90x90  matrix  of  the  minimum  genetic  distance  among  breeding 
groups  was  estimated  at  year  2005  in  each  simulation,  to  reflect  the  accumulated 
influence  gene  flow  and  genetic  drift  had  on  the  genetic  differences  among  groups. 

We  used  a  similar  approach  to  estimate  the  observed  minimum  genetic  distance  among 
territories.  The  MCBCL  monitoring  data  were  used  to  compile  a  pedigree  of  all  breeding 
events  from  1986  to  2005.  A  gene  drop  was  then  applied  to  this  pedigree  assuming  an 
infinite  alleles  model.  At  year  2005  the  minimum  genetic  distance  was  estimated.  This 
process  was  repeated  10,000  times  and  the  estimated  observed  minimum  genetic  distance 
represents  an  average  of  this  process. 

A  Mantel  test  was  used  to  test  for  the  degree  of  similarity  between  the  observed 
minimum  genetic  distance  and  minimum  genetic  distance  generated  in  the  simulation 
(Bailey  and  Gatrell  1995).  An  average  Mantel  correlation  coefficient  from  the  50 
replicates  was  used  to  characterize  the  ability  of  each  parameterization  to  reproduce 
genetic  patterns  observed  on  MCBCL,  MrD0bs. 

To  derive  the  pattern  fulfillment  criterion,  we  calculated  a  distribution  of  Mantel  test 
statistics  by  randomly  exchanging  minimum  genetic  distance  values  among  territories 
(i.e.,  random  permutation)  and  correlating  this  random  pattern  with  the  observed  pattern, 
Mr  Rand-  This  was  performed  10,000  times  as  explained  above  for  Group  Size  and 
Connectivity.  The  pattern  fulfillment  criterion  for  connectivity  specified  that 
parameterizations  must  produce  Mr°obs  greater  than  0.069  (p<0.0001). 

1.3. 5. 6  Filter  Optimization.  For  model  parameterizations  capable  of  fulfilling  all 
patterns,  we  investigated  the  sensitivity  of  pattern  fulfillment  indices  to  the  remaining 
parameter  uncertainty.  This  was  accomplished  by  performing  a  Spearman  rank 
correlation  between  the  remaining  range  of  possible  parameter  values  and  each  of  the  five 
filters.  Large  correlation  coefficients  indicate  that  model  fit  could  be  improved  by  further 
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restricting  pattern  fulfillment  criteria.  Filter  optimization  also  indicates  which  patterns 
are  most  useful  for  reducing  prediction  uncertainty  and  can  be  used  to  guide  further  data 
collection. 


1.4  Results  and  Discussion 

Parameterizations  generating  patterns  unable  to  meet  all  five  pattern  fulfillment  criteria  were 
excluded  from  further  analysis.  We  found  that  zero  parameterizations  produced  patterns  meeting 
all  five  criteria.  Examining  the  distributions  of  pattern  fulfillment  indices  indicated  that  only  8 
parameterizations  produced  Emf  loater  less  than  or  equal  to  2.9.  Of  these  8  parameterizations,  none 
were  able  to  meet  the  pattern  fulfillment  criteria  for  the  remaining  filters.  Recall  that  very  few 
male  floaters  are  observed  in  the  landscape  (Fig.  1.3)  and  that  we  expect  a  large  amount  of 
observational  error  associated  with  this  monitoring  data.  We  decided  to  allow  for  a  larger 
amount  of  error  for  this  pattern  fulfillment  criteria  by  specifying  that  Emlloatcr  must  be  less  than  or 
equal  to  90%  0f  p;mFloater  observed  in  the  100,000  parameterizations,  Eml  loater<=  7.32. 

After  applying  the  adjusted  pattern  fulfillment  criteria  for  male  floaters,  we  found  that  only  719 
parameterizations  produced  patterns  meeting  all  five  criteria.  Therefore,  filtering  removed 
99.3%  of  all  parameter  sets  tested — a  significant  reduction  in  uncertainty.  However,  the 
parameter  state  space,  or  range  of  parameter  values,  for  the  2 1  uncertain  parameters  remaining 
after  filtering  often  equaled  the  initial  parameter  state  space  included  in  the  simulations. 
Therefore,  the  percent  of  parameter  state  space  removed  was  often  zero  (Table  1.4).  The  biggest 
exception  was  observed  for  probability  of  male  floater  survival,  indicating  the  35.9%  of  the 
parameter  state  space  has  been  removed.  Table  1.4  also  indicates  that  -50%  of  parameter  state 
space  for  d45°,  d9o°,  and  di35°  was  removed.  This  filtering  actually  reflects  the  autocorrelation  in 
the  probability  of  turning  described  above  in  the  Habitat-based  Correlated  Random  Walk  model. 

These  results  indicate  that  the  Letcher  model  is  insensitive  to  the  dispersal  rules  included  in  the 
model.  In  other  words,  the  demographic  parameters  and  behavioral  rules  that  were  not  treated  as 
uncertain  are  capable  of  approximating  patterns  observed  on  MCBCL.  These  results  may  also 
reflect  the  clumped  spatial  distribution  of  RCW  territories  on  MCBCL.  In  other  words, 
individuals  do  not  have  to  travel  far  or  interact  much  with  landscape  structure  in  order  to 
facilitate  the  exchange  of  individuals  among  groups  required  to  reproduce  patterns. 
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Table  1.4.  Results  of  filtering  for  the  21  uncertain  parameters  using  the  five  pattern  fulfillment 
criteria.  N=719  (number  of  parameterizations  remaining  after  filtering). 

Standard  Percent1 

Uncertain  Parameter  Mean  Deviation  filtered 


fFS 

Seasonal  probability  of  survival  -  Female  Floater 

0.840 

0.034 

0 

mFS 

Seasonal  probability  of  survival  -  Male  Floater 

0.824 

0.023 

35.9 

mHP 

Perceptual  distance  of  Male  Helpers 

8.51 

12.7 

0 

mFP 

Perceptual  distance  of  Male  Floaters 

38.7 

13.5 

13.6 

fHP 

Perceptual  distance  of  Female  Helpers 

24.1 

17.1 

0 

fFP 

Perceptual  distance  of  Female  Floaters 

25.1 

12.4 

13.6 

mSp 

Male  floater  dispersal  speed 

43.8 

21.0 

0 

fSp 

Female  floater  and  Male  natal  dispersal  speed 

56.6 

27.2 

0 

do° 

Probability  of  traveling  in  same  direction  that  was 
chosen  in  last  time  step 

0.532 

0.304 

0 

d45° 

Probability  of  turning  45°  away  from  direction 
traveled  in  previous  time  step 

0.122 

0.114 

50.0 

d90° 

Probability  of  turning  90°  away  from  direction 
traveled  in  previous  time  step 

0.058 

0.071 

55.0 

dl35° 

Probability  of  turning  135°  away  from  direction 
traveled  in  previous  time  step 

0.054 

0.072 

55.8 

ml 

Probability  of  entering  predominantly  pine  forest 

0.506 

0.306 

0 

m2 

Probability  of  entering  mixed  pine  and  hardwood 
forest 

0.501 

0.302 

0 

m3 

Probability  of  entering  bottomland  hardwood  / 
upland  hardwood  forest 

0.516 

0.301 

0 

m4 

Probability  of  entering  forest  plantation  under  10 
years 

0.520 

0.302 

0 

m5 

Probability  of  entering  bare  ground  /  cropland  / 
grasses  /  marshland  /  scrub,  shrub 

0.488 

0.298 

0 

m6 

Probability  of  entering  river,  stream,  lake 

0.526 

0.305 

0 

m7 

Probability  of  entering  business  or  commercial  area 
/  residential  area 

0.505 

0.306 

0 

tn8 

Probability  of  entering  off  site  areas 

0.497 

0.290 

0 

B 

Strength  of  preference  for  habitat-based  or 
directional  dispersal  rules  (e.g.,  B  =  1,  only  land 
cover  type  matters) 

0.475 

0.312 

0 

'Percent  filtered  =  1  -  [range  of  values  from  719  /  range  of  values  from  100,000  parameterizations] 
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Filter  Optimization.  Spearman  rank  correlations  between  parameter  values  remaining  after 
filtering  and  pattern  fulfillment  indices  indicated  that  all  filters  are  sensitive  to  the  remaining 
uncertainty  in  parameter  values  (i.e.,  Table  1.5,  p<0.05  observed  for  at  least  one  parameter  in  all 
indices).  Therefore,  by  further  restricting  the  pattern  fulfillment  criteria  we  should  be  able  to 
reduce  a  considerable  amount  of  uncertainty  in  these  parameters.  Spearman  rho  values  indicate 
that  Mantel’s  r  for  female  and  male  connectivity,  Error  PBG,  and  Pearson’s  r  for  Group  Size  are 
strong  filters,  indicating  that  continuation  of  the  current  monitoring  program  will  help  reduce 
uncertainty  regarding  dispersal  behaviors. 

The  largest  correlations  were  observed  for  female  floater  perceptual  distance.  Female  floater 
perceptual  distance  was  negatively  correlated  with  EPBG  (rho  =  -0.736).  Therefore,  the  longer  a 
female  floater’s  perceptual  distance,  the  error  in  predicting  PBG  will  be  smaller.  Female  floater 
perceptual  distance  was  negatively  correlated  with  Mantel’s  r  for  female  connectivity  (-0.768). 
Therefore,  the  longer  a  female  floater’s  perceptual  distance,  the  correlation  between  observed 
and  expected  female  connectivity  will  be  smaller  (Fig  1.4).  Further,  the  relationship  between 
female  floater  perceptual  distance  and  Mantel’s  r  for  minimum  genetic  distance  shows  a  positive 
relationship  (0.549).  Therefore,  the  correlation  between  observed  and  expected  minimum 
genetic  distance  among  clusters  is  greater  when  female  floater  perceptual  distance  is  larger.  By 
using  multiple  filters  we  will  likely  find  an  intermediate  value  for  female  floater  perceptual 
distance  that  would  be  otherwise  missed  if  only  one  filter  were  used. 


Female  Floater  Perceptual  Distance 


Figure  1,4.  Scatter  plot  of  parameter  values 
for  female  floater  perceptual  distance  and 
corresponding  pattern  fulfillment  criteria 
(N=719).  Significant  Spearman’s  rho  (Table 
1.5)  indicates  that  restricting  pattern 
fulfillment  criteria  will  help  further  reduce 
parameter  uncertainty. 
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Table  1.5.  Spearman  Rank  correlation  coefficient  (rho)  describing  the  strength  of  the  relationship  between  uncertain  parameter  values 
and  statistical  indices  of  model  fit  to  observed  patterns  on  MCBCL.  N=719. 

r  Mr  Mr 


Parameter 

Error 

PBG 

P- 

value 

Error 

mFloat 

P- 

value 

Error 

fFloat 

P- 

value 

Group 

Size 

P- 

value 

Male 

Conn 

P- 

value 

Female 

Conn 

P- 

value 

Mr  Dij 

P- 

value 

fFS 

0.058 

0.123 

-0.030 

0.419 

0.324 

<0.001 

0.018 

0.622 

0.070 

0.061 

-0.070 

0.061 

0.014 

0.710 

mFS 

0.155 

<0.001 

0.065 

0.080 

-0.098 

0.008 

-0.038 

0.308 

0.026 

0.494 

0.246 

<0.001 

-0.028 

0.453 

mHP 

0.446 

<0.001 

0.292 

<0.001 

-0.127 

0.001 

-0.145 

<0.001 

-0.661 

<0.001 

0.381 

<0.001 

-0.353 

<0.001 

mFP 

-0.133 

<0.001 

-0.075 

0.043 

0.134 

<0.001 

0.653 

<0.001 

-0.554 

<0.001 

0.235 

<0.001 

0.135 

<0.001 

fHP 

-0.232 

<0.001 

0.033 

0.378 

-0.073 

0.051 

-0.001 

0.969 

0.049 

0.193 

-0.143 

<0.001 

0.003 

0.937 

fFP 

-0.738 

<0.001 

0.167 

<0.001 

0.414 

<0.001 

0.384 

<0.001 

0.609 

<0.001 

-0.768 

<0.001 

0.549 

<0.001 

mSp 

0.023 

0.535 

-0.037 

0.318 

-0.097 

0.009 

-0.009 

0.818 

-0.036 

0.339 

0.022 

0.557 

0.002 

0.959 

fSp 

-0.047 

0.204 

-0.030 

0.415 

0.000 

0.994 

0.036 

0.338 

-0.038 

0.310 

-0.174 

<0.001 

0.065 

0.080 

d0” 

-0.058 

0.123 

-0.048 

0.201 

-0.002 

0.962 

0.081 

0.031 

-0.014 

0.700 

-0.127 

0.001 

0.090 

0.016 

ml 

0.020 

0.597 

0.023 

0.545 

-0.059 

0.114 

-0.012 

0.745 

0.008 

0.829 

0.056 

0.137 

0.019 

0.617 

m2 

-0.030 

0.417 

-0.047 

0.209 

0.026 

0.480 

0.039 

0.292 

0.041 

0.270 

0.018 

0.629 

0.066 

0.078 

m3 

0.033 

0.381 

0.027 

0.462 

0.037 

0.327 

-0.027 

0.465 

0.002 

0.948 

0.024 

0.524 

-0.068 

0.068 

m4 

0.031 

0.406 

-0.007 

0.847 

-0.020 

0.593 

-0.034 

0.359 

0.003 

0.942 

0.050 

0.177 

0.049 

0.189 

m5 

0.001 

0.989 

0.005 

0.892 

0.036 

0.339 

0.009 

0.809 

0.015 

0.690 

-0.004 

0.921 

0.043 

0.253 

m6 

-0.026 

0.489 

-0.033 

0.373 

0.030 

0.420 

0.064 

0.084 

-0.081 

0.030 

0.066 

0.078 

-0.010 

0.793 

m7 

0.041 

0.277 

-0.039 

0.295 

-0.016 

0.665 

-0.035 

0.355 

-0.028 

0.448 

0.038 

0.314 

-0.024 

0.528 

m8 

0.004 

0.915 

0.015 

0.684 

0.063 

0.092 

0.051 

0.169 

-0.055 

0.141 

-0.001 

0.973 

0.016 

0.667 

B 

-0.070 

0.061 

0.088 

0.018 

0.062 

0.097 

0.051 

0.169 

0.086 

0.021 

0.129 

0.001 

0.066 

0.076 

Bold  values  =  p<0.05 
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Spearman  rho  values  indicate  that  parameter  values  of  the  Habitat-based  Correlated  Random 
Walk  tended  to  be  insensitive  to  the  five  pattern  fulfillment  indices.  We  did  observe  negative 
correlations  between  the  probability  of  not  turning  (do  )  and  model  fit  for  female  connectivity 
(p=0.001).  This  suggests  that  flexibility  in  directional  choices  is  important  for  approximating 
connectivity  values  observed  in  the  landscape.  In  contrast,  a  small  positive  correlation  was 
observed  between  Mantel’s  r  for  minimum  genetic  distance  and  the  probability  of  not  turning 
(rho  =  0.09,  p  =  0.016).  Therefore,  straighter  dispersal,  which  would  allow  individuals  to  move 
farther,  does  improve  the  fit  between  observed  and  expected  genetic  distance.  Among  the 
preferences  for  land  cover  types,  only  m6,  probability  of  crossing  water,  shows  a  correlation  with 
Mantel’s  r  for  male  connectivity  (rho  =  -0.081,  p  =  0.03).  The  strength  of  preference  for  habitat- 
based  dispersal  rules  (B)  was  positively  correlated  with  model  fit  for  male  and  female 
connectivity.  In  contrast,  EmF  oat  was  positively  correlated  with  preference  for  habitat-based 
dispersal  rules.  Therefore,  when  B  is  higher  we  observed  a  poorer  fit  between  observed  and 
expected  number  of  male  floaters  in  the  landscape.  This  correlation  may  indicate  that  habitat 
based  dispersal  rules  lead  to  retaining  more  floaters  in  the  landscape  than  straight  dispersal, 
which  leads  to  greater  emigration  away  from  MCBCL. 

We  investigated  how  further  constraining  pattern  fulfillment  criteria  may  reduce  parameter 
uncertainty.  We  derived  4  filter  combinations  by  increasing  strength,  or  decreasing  filter 
porosity,  of  the  pattern  fulfillment  criteria.  Derivation  filter  combination  A  was  explained  above 
using  random  permutation  testing  and  led  to  the  reduced  719  parameterizations  (Table  1.6). 

Using  the  distribution  of  pattern  fulfillment  indices  that  characterize  the  719  parameterizations, 
we  specified  that  the  pattern  fulfillment  criteria  must  be  less  than  the  lower  50th  percentile  for 
Epbg  and  greater  than  the  50th  percentile  of  correlation  coefficients  for  group  size,  connectivity, 
and  minimum  genetic  distance.  Twenty  parameterizations  were  able  to  satisfy  all  criteria.  This 
process  was  repeated  by  adjusting  the  criteria  to  be  lower  than  the  40th  and  greater  than  the  upper 
60th  percentile,  leaving  8  parameterizations  capable  of  meeting  all  five  criteria.  Subsequently 
adjusting  for  the  lower  30th  and  upper  70th  percentile,  we  found  only  one  parameterization  met  all 
five  criteria  (Table  1.6). 
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Table  1.6.  Filter  combinations  used  to  reduce  parameter  uncertainty 


Pattern 

Filter 

Filter 

Filter 

Filter 

Fulfillment 

Combination 

Combination 

Combination 

Combination 

Index 

A 

B 

C 

D 

Error  PBG 

15.2 

9.58 

8.26 

7.02 

'Error  Male 
Floater 

7.32 

7.32 

7.32 

7.32 

1  Error  Female 
Floater 

Pearson's  r 

7.95 

7.95 

7.95 

7.95 

Group  Size 
Mantel's  r  Male 

0.558 

0.659 

0.662 

0.667 

Connectivity 
Mantel's  r 

0.116 

0.254 

0.257 

0.261 

Female 
Connectivity 
Mantel's  r 

Genetic  Distance 

0.065 

0.097 

0.1 

0.103 

(Dij) 

Number  of 

0.069 

0.173 

0.178 

0.182 

parameterizations 

remaining 

719 

20 

8 

1 

Filter  Combination  A  =  derived  from  random  permutation  testing 
Filter  Combination  B  =  represents  lower  and  upper  50th  percentiles 
Filter  Combination  C  =  represents  lower  40th  and  upper  60th  percentiles 
Filter  Combination  D  =  represents  lower  30th  and  upper  70th  percentiles 

'E1'1""  was  not  reduced  for  filter  combinations  B,  C  and  D  because  in  Filter  A  EmFloat  equaled  the  lower  10th 
percentile  and  Effloat  is  less  than  the  30th  percentile. 

Figure  1.5  illustrates  how  the  parameter  estimates  and  uncertainty  change  with  increasing  the 
strength  of  the  pattern  fulfillment  criteria.  A  considerable  amount  of  uncertainty  was  removed 
regarding  floater  survival.  The  filter  most  directly  related  to  floater  survival  (i.e.,  EFloat)  was  not 
changed  in  the  successive  filter  combinations  (Table  1.6),  therefore  the  mean  value  does  not  shift 
much.  In  contrast,  perceptual  distance  for  males  did  change  considerably  with  filter 
combinations.  POM  predicts  that  male  floaters  have  a  much  greater  perceptual  distance  than 
male  helpers.  For  females,  POM  predicts  that  perceptual  distance  for  helpers  and  floaters 
converge  to  the  same  value.  These  results  suggest  that  male  helpers  and  floaters  engage  in  very 
different  behaviors  while  female  helpers  and  floaters  have  more  similar  behaviors.  These  results 
seem  to  agree  with  observations  of  Jeffrey  Walter’s  field  crews  conducting  the  radio  telemetry 
study  of  RCWs  on  Fort  Bragg  (Jeffrey  Walters,  personal  communication). 

Uncertainty  regarding  male  dispersal  speed  was  not  reduced  in  the  successive  filter  combinations 
and  a  large  shift  was  observed  with  filter  combination  D.  However,  because  male  floater 
perceptual  distance  was  high  in  the  simulations  passing  through  these  filters,  male  dispersal 
speed  is  likely  not  a  critical  parameter  required  to  approximate  patterns  observed  in  the  field.  In 
contrast,  uncertainty  regarding  female  dispersal  speed  was  reduced  and  the  mean  shifted  lower 
with  successive  filter  combinations. 
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Uncertainty  regarding  the  parameter  values  for  the  HCRW  was  not  reduced  with  successive  filter 
combinations,  and  large  shifts  in  parameter  estimates  are  often  observed  with  Filter  Combination 
D.  For  example  for  m4,  probability  of  entering  a  cell  containing  young  pine,  the  mean  shifts 
upward  until  inclusion  of  filter  combination  D.  Similarly,  for  m5,  probability  of  entering  a  non- 
forested  /  undeveloped  cell,  the  mean  shifts  lower  until  inclusion  of  filter  combination  D. 
Therefore,  we  can  not  place  much  confidence  in  the  parameter  values  derived  from  filter 
combination  D  for  the  HCRW. 
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Figure  1.5.  Change  in  parameter  estimates  for  the  uncertain  dispersal  parameters  owing  to  filter 
optimization.  Error  bars  represent  one  standard  deviation.  Filter  Comb  A,  n=719.  Filter  Comb  B, 
n=20.  Filter  Comb  C,  n=8.  Filter  Comb  D,  n=l  (Table  1.6). 
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1.5  Conclusions 

POM  did  reduce  a  considerable  amount  of  uncertainty  in  the  RCW  SEPM,  especially  for 
parameters  relating  to  floater  survival  and  perceptual  distance.  However,  POM  was  unable  to 
reduce  uncertainty  regarding  the  HCRW.  One  reason  may  be  that  the  clumped  distribution  of 
habitat  allows  the  SEPM  to  approximate  patterns  observed  in  the  field  regardless  of  the  dispersal 
rules  employed.  Another  reason  may  be  that  there  is  structural  uncertainty  in  the  model. 
Structural  uncertainty  refers  to  how  the  submodels  simulate  biological  processes  and  how  the 
submodels  are  interconnected. 

While  the  HCRW  has  been  successfully  applied  to  other  species  (Kramer- Schadt  et  al.  2004), 
assumptions  used  in  other  submodels  originally  included  in  the  Letcher  model  may  prevent  its 
usefulness  for  RCWs.  For  example,  the  Letcher  competition  submodel  indicates  that  all  floaters 
and  helpers  within  3  km  (or  30  grid  cells  in  our  simulations)  compete  for  a  breeding  vacancy. 
Recall  that  this  distance  reflects  the  bird’s  perceptual  distance.  Within  a  bird’s  perceptual 
distance,  the  intervening  land  cover  types  are  assumed  to  have  no  influence  on  movement, 
keeping  with  the  basic  structural  assumptions  of  the  Letcher  model. 

We  tested  a  bird’s  perceptual  distance  with  POM,  by  allowing  perceptual  distance  to  vary 
independently  for  helpers  and  floaters  and  detennining  the  distances  that  allow  the  simulation 
model  to  approximate  multiple  patterns  observed  in  the  field.  For  males,  we  found  that  helpers 
only  compete  for  vacancies  very  close  to  their  location  (i.e.,  within  300  m)  and  floaters  can 
compete  for  vacancies  very  far  away  from  their  current  location  (i.e.,  within  4.5  km)  (Figure 
1 .5).  In  contrast,  for  females  we  found  that  both  helpers  and  floaters  can  compete  for  vacancies 
very  far  from  their  current  location  (~2.4  km). 

Perceptual  distance  is  the  main  parameter  in  the  competition  submodel.  The  competition 
submodel  seems  to  be  the  dominant  component  of  the  RCW  SEPM  that  determines  patterns 
observed  on  MCBCL.  The  HCRW  simulates  “searching  behaviors”  of  floaters  when  no 
breeding  vacancy  is  detected  within  their  perceptual  distance.  The  interaction  between  searching 
behaviors  and  competition  requires  further  study.  Searching  behaviors  are  likely  to  contribute 
more  to  population  viability  when  habitats  become  more  fragmented.  The  habitat  trading 
scenarios  examined  in  the  subsequent  section  may  increase  habitat  fragmentation.  Therefore,  we 
chose  to  evaluate  these  scenarios  with  the  8  parameterizations  that  passed  through  filter 
combination  C.  This  approach  allows  us  to  incorporate  the  remaining  uncertainty  in  the  HCRW 
by  including  multiple  models  of  the  influence  of  land  cover  on  dispersal  behaviors. 
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2.  Phase  2.  Habitat  Trading  for  the  Red-Cockaded  Woodpecker  under 
Uncertainty  Regarding  Dispersal  Behaviors 


2.1  Objective 

The  objective  of  our  SEED  project  was  to  derive  a  method  of  making  decisions  regarding  habitat 
trades  while  incorporating  uncertainty  regarding  how  landscape  pattern  affects  dispersal 
behaviors,  and  ultimately  landscape  connectivity.  First,  we  described  a  landscape-scale  tradable 
credit  system  designed  to  prevent  habitat  trades  that  exacerbate  habitat  fragmentation.  Second, 
we  applied  the  tradable  credit  system  to  hypothetical  trading  scenarios  for  MCBCL.  The  trading 
scenarios  are  evaluated  using  the  Letcher  model  currently  used  for  RCW  management  by 
USFWS  and  DoD.  We  also  used  the  eight  Habitat-based  Correlated  Random  Walk  models 
derived  from  POM  in  Phase  1  to  evaluate  habitat  trades.  Third,  we  applied  the  principles  of 
decision  analysis  to  weigh  alternative  management  decisions  in  the  face  of  uncertainty.  By 
integrating  Pattern-Oriented  Modeling  (Phase  1)  with  a  landscape-scale  tradable  credit  system, 
we  described  how  one  can  determine  when  uncertainty  is  sufficiently  large  that  collecting  more 
data  and  reducing  uncertainty  is  a  better  choice  than  trading  habitat.  Conversely,  the  approach 
can  be  used  to  determine  when  enough  information  has  been  collected  and  uncertainty  regarding 
dispersal  behaviors  is  sufficiently  small  to  allow  habitat  trading. 

2.2  Background 

Marine  Corp  Base  Camp  Lejeune  (MCBCL)  is  one  of  three  management  units  for  the  Coastal 
North  Carolina  Population  of  Red-Cockaded  Woodpeckers  (RCWs).  The  other  two  include 
Croatan  National  Forest  and  Holly  Shelter  Game  Lands.  All  three  management  units  are  located 
within  the  Onslow  Bight  Landscape,  which  stretches  from  Cape  Lookout  to  Cape  Fear  (Fig  2.1). 
The  Onslow  Bight  Partnership  consists  of  federal  and  non-federal  natural  resource  managers 
dedicated  to  preserving  the  biodiversity  of  the  region.  At  recovery,  the  Coastal  North  Carolina 
population  will  include  350  Potential  Breeding  Groups  (PBGs),  which  includes  an  adult  male 
and  female  with  or  without  helpers  and  fledglings  (USFWS  2003a).  To  contribute  to  recovery, 
the  USFWS  has  specified  that  MCBCL  maintain  a  population  growth  rate  of  5%  per  year  and 
achieve  173  active  clusters  by  2025.  Active  clusters  are  defined  as  those  which  contain  a  solitary 
individual,  may  be  captured  by  a  neighboring  group,  or  contain  a  PBG  (USFWS  2003a).  When 
MCBCL  reaches  the  173  active  cluster  goal,  training  restrictions  due  to  RCWs  will  be  eliminated 
or  relieved  (MCBCL  2007).  Establishing  new  clusters  outside  of  MCBCL,  Croatan,  and/or 
Holly  Shelter  within  the  Onslow  Bight  Landscape,  through  DoD’s  Encroachment  Partnering 
program,  will  speed  when  regulatory  relief  is  achieved  (TNC  2006;  Ralph  Costa,  USFWS, 
personal  communication,  2006). 

2.2.1  Habitat  Trading 

A  quantitative  method  for  detennining  incremental  changes  in  regulatory  relief  thanks  to 
Encroachment  Partnering,  prior  to  achieving  population  recovery  goals,  would  benefit  many 
DoD  facilities.  Due  to  the  Base  Realignment  and  Closure  process,  training  on  many  DoD 
facilities  has  intensified.  Below  we  discuss  two  important  implications  of  using  the 
Encroachment  Partnering  program  to  achieve  regulatory  relief  for  DoD.  The  first  is  that  the 
Endangered  Species  Act  (USFWS  1973)  applies  different  standards  to  non-federal  and  federal 
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Figure  2.1.  The  Onslow  Bight  Landscape  in  southern  coastal  North  Carolina. 
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landowners.  The  second  is  that  achieving  population  recovery  through  partnerships  among 
multiple  landowners  has  the  capacity  to  drastically  change  where  habitats  are  located  for  listed 
species.  The  connectivity  of  habitats  plays  an  important  role  in  the  recovery  and  long-term 
viability  of  many  species  (Hanski  and  Gaggiotti  2004). 

The  different  expectations  for  federal  and  non-federal  landowners  are  found  in  Sections  7,  9,  and 
10  of  the  ESA  (USFWS  1973).  Non-federal  landowners  are  prevented  from  taking  endangered 
species,  which  is  defined  as  an  action  that  harms,  harasses,  or  kills  a  listed  species  (Section  9; 
USFWS  1973).  Adverse  modification  to  habitat  that  results  in  death  or  injury  is  included  in  the 
take  provision.  However,  non-federal  landowners  are  not  required  to  manage  their  lands  for 
listed  species,  and  loss  of  natural  disturbance  regimes  (e.g.,  fire)  is  often  enough  to  decrease  the 
probability  of  habitat  occupancy  (James  et  al.  2001).  Under  Section  10,  non-federal  landowners 
may  take  individuals  or  remove  habitat  if  a  Habitat  Conservation  Plan  is  implemented  to  offset 
the  take  elsewhere  in  the  landscape  (USFWS  1973).  This  has  led  to  the  adoption  of 
“conservation  banking,”  or  a  tradable  credit  system  for  endangered  species  habitat  (USFWS 
2003b).  Tradable  credit  systems  represent  the  trading  of  access  rights  to  endangered  species 
habitat.  Note  that  federal  landowners  are  prevented  from  providing  habitat  to  offset  the  loss  of 
habitat  on  non-federal  lands  (USFWS  2003b).  This  reflects  the  Section  7  requirement  that 
federal  landowners  must  pro-actively  contribute  to  the  recovery  of  listed  species  (USFWS  1973). 
Section  7  further  prevents  federal  landowners  from  any  action  that  would  appreciably  reduce  the 
likelihood  of  survival  and  recovery,  or  jeopardize  the  continued  existence  of  the  species 
(Stanford  Environmental  Law  Society  2001). 

The  RCW  Recovery  Plan  discusses  recovery  goals  and  conditions  required  for  offsetting  a  take 
in  terms  of  PBGs  and  active  clusters  (USFWS  2003a),  supporting  the  broader  regulatory  goal  of 
preventing  the  net  loss  of  habitat  area  for  endangered  species  (USFWS  1973).  However, 
changing  the  connectivity  of  RCW  clusters  can  also  affect  RCW  population  viability  (Cox  et  al. 
2001;  Letcher  et  al  1998).  Landscape  connectivity  can  be  defined  as  “the  degree  to  which  the 
landscape  facilitates  or  impedes  movement  among  resource  patches”  (Taylor  et  al.  1993).  While 
there  are  no  regulatory  criteria  for  managing  changes  in  habitat  connectivity  in  the  ESA,  the 
RCW  Recovery  Plan  does  consider  this  problem  under  the  Private  Lands  Conservation  Strategy 
(USFWS  2003a).  RCW  Recovery  Plan  specifies  that  only  isolated  clusters,  defined  as  those 
sharing  few  migrants  with  other  clusters,  can  be  traded  for  new  clusters  that  are  more  connected. 
Further,  the  location  where  new  clusters  are  established  must  have  capacity  for  at  least  10  highly 
connected  RCW  groups.  However,  defining  levels  of  connectivity  at  which  habitat  trading  is  a 
viable  management  option  is  an  important  area  of  research  (USFWS  2003a).  To  achieve  military 
readiness,  the  DoD  will  often  wish  to  change  the  allocation  of  habitat  for  endangered  species 
prior  to  reaching  scientific  consensus  regarding  dispersal  behaviors. 

2.2.2  Considering  Landscape  Processes  when  Trading 

Parcels  will  differ  in  their  contribution  to  total  (meta)population  growth,  rates  of  migration,  and 
maintenance  of  genetic  variance,  based  not  only  on  the  size  and  quality  of  the  habitat  patch,  but 
also  the  landscape  context  of  a  parcel  (i.e.,  its  location  relative  to  all  other  land  uses).  For 
example,  even  if  habitat  patches  of  equivalent  size,  shape,  and  quality  are  traded,  habitat  trading 
will  likely  change  the  rates  of  migration  among  other  patches  within  the  region.  The  change  in 
regional,  landscape  structure  may  subsequently  change  rates  of  local  extinction  and  erosion  of 
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genetic  variance  in  subdivided  populations  such  as  RCWs  (Bruggeman  2005).  A  measure  of 
conservation  value  that  captures  the  interaction  between  local  and  regional  biological  processes 
is  needed  to  ensure  that  an  ecologically  equivalent  trade  in  habitat  is  made  (Bruggeman  et  al. 
2005).  A  tradable  credit  system  that  accounts  for  differences  in  conservation  value  among 
habitat  patches  is  more  likely  to  achieve  sustainable  populations  at  the  least  cost. 

We  argue  that  that  population  genetic  theory  provides  a  robust  approach  for  avoiding  trades  that 
increase  habitat  fragmentation  and,  subsequently,  decrease  population  viability.  Subdivided 
populations  must  achieve  a  balance  between  growth  within  each  breeding  group  and  migration 
between  breeding  groups  to  prevent  local  extinction,  inbreeding  and  outbreeding  depression,  and 
to  maintain  genetic  variance  needed  for  adaptation  through  natural  selection  at  both  individual 
and  group  levels  (Harrison  and  Hastings  1996;  Mills  and  Allendorf  1996).  Neutral  genetic 
variance  is  useful  for  determining  how  habitat  loss  or  restoration  affects  gene  flow,  genetic  drift, 
and  inbreeding  (Hedrick  2001).  Field  studies  have  demonstrated  that  mammalian  populations 
have  evolved  behaviors  to  simultaneously  prevent  loss  of  genetic  variance  within  while 
maintaining  genetic  variance  among  local  populations  (Dobson  et  al.  1997;  Storz  1999;  Coltman 
et  al.  2003).  Using  a  computer  simulation,  Bruggeman  (2005)  showed  that  the  complex  mating 
system  of  RCWs  (as  described  by  Walters  1988)  is  very  efficient  at  maintaining  genetic  variance 
within  and  among  breeding  groups  when  a  continuous  distribution  of  habitat  exists,  but  habitat 
fragmentation  increased  rates  of  inbreeding  and  local  extinction.  In  a  population  of  220  RCW 
groups,  Daniels  and  Walters  (2000)  observed  that  short  dispersal  distances  can  result  in  the 
mating  between  closely  related  individuals,  causing  inbreeding  depression,  measured  as  a 
reduction  in  recruitment. 

The  spatial  distribution  of  genetic  variance  is  a  reflection  of  the  rates  of  growth  within  local 
populations,  migration,  and  local  extinction.  The  spatial  apportionment  of  genetic  variance 
observed  in  a  landscape  with  an  allocation  of  habitat  observed  prior  to  habitat  loss  and 
fragmentation,  a  baseline  condition,  is  indicative  of  rates  of  migration  and  population  growth  that 
maintain  adaptive  genetic  variance  and  prevent  inbreeding  and  outbreeding  depression. 

Therefore,  we  can  quantify  changes  in  fragmentation  by  comparing  the  spatial  apportionment  of 
genetic  variance  before  and  after  the  trade  to  that  observed  in  this  baseline  condition  (Bruggeman 
et  al.  2005).  If  trading  leads  to  greater  genetic  variance  among  local  populations,  compared  to 
this  baseline  condition,  then  rates  of  migration,  recolonization,  and/or  local  population  growth 
have  decreased  and  genetic  drift,  inbreeding  depression,  and  local  extinction  may  threaten 
viability.  If  trading  leads  to  greater  genetic  variance  within  local  populations  then  outbreeding 
depression  and  lost  opportunities  for  local  adaptations  may  threaten  viability. 

2.2.3  Landscape  Equivalency  Analysis  (LEA) 

By  integrating  population  genetic  theory  and  natural  resource  economics,  we  have  devised  a 
tradable  credit  system  to  incorporate  the  influence  of  changing  habitat  area  and  connectivity  on 
the  viability  of  endangered  species  (Bruggeman  et  al.  2005).  We  developed  Landscape 
Equivalency  Analysis  (LEA)  as  a  generally  applicable,  landscape-scale  accounting  system  to 
calculate  the  conservation  value  of  tradable  credits  for  endangered  species  habitat  to  minimize 
and  reverse  the  negative  effects  of  habitat  fragmentation  (Bruggeman  et  al.  2005).  Endangered 
species  management  may  benefit  from  incorporating  habitat  defragmentation  into  recovery 
goals.  Therefore,  LEA  adds  a  “spatially-explicit”  recovery  objective  as  the  allocation  of  habitat 
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yielding  the  spatial  apportionment  of  genetic  variance  (e.g.,  as  would  be  determined  using 
neutral  genetic  markers)  observed  prior  to  habitat  loss  and  fragmentation  (Meffe  1996; 
Bruggeman  et  al.  2005).  Given  existing  loss  of  habitat  in  most  landscapes,  this  goal  will  often  be 
unachievable,  but  it  does  serve  as  an  objective  criterion  for  defragmenting  endangered  species 
habitat  to  protect  evolutionary  processes  at  the  landscape  scale. 

LEA  is  an  extension  of  resource-based  compensation  applied  to  a  landscape-scale.  When 
environmental  regulations  stipulate  in-kind  replacement  of  ecological  resources,  as  in  the  ESA, 
compensation  must  be  made  using  the  same  type  of  resources  that  were  lost  (Jones  and  Pease 
1997;  Zafonte  and  Hampton  2007).  Resource-based  compensation  determines  the  amount  of 
restoration  required  to  equate  an  individual’s  well-being  before  loss  of  an  ecological  resource 
with  their  well-being  after  habitat  destruction.  A  “service-to-service”  approach  is  used  wherein 
adequate  compensation  is  made  to  the  public  if  the  habitats  restored  provide  equivalent  types  and 
levels  of  ecological  services  as  the  habitats  lost.  Ecological  functions  can  be  treated  as  goods 
and  services  when  a  direct  or  indirect  benefit  to  humans  can  be  demonstrated  (deGroot  et  al. 
2002).  The  direct  human  benefits  of  protecting  endangered  species  include  use  value  (e.g., 
seeing  the  species),  option  value  (e.g.,  possibility  that  genetic  variance  provided  by  the  species 
may  contribute  to  medical  or  agricultural  advances),  existence  value  (i.e.,  knowing  the  species 
exists),  and  bequest  value  (i.e.,  knowing  the  species  will  be  present  for  future  generations) 
(Loomis  and  White  1996).  Genetic  variance  of  an  endangered  population  also  provides  indirect 
benefits  for  humans.  Adaptive  genetic  variance  is  required  for  population  persistence  (Fisher 
1930).  As  previously  discussed,  neutral  genetic  variance  provides  useful  information  regarding 
how  populations  respond  to  changes  in  habitat  area  and  connectivity.  Thus,  endangered  species 
habitat  is  a  form  of  natural  capital,  defined  as  a  stock  of  resources  providing  useful  services 
(deGroot  et  al.  2000).  In  this  analysis  we  assume  that  the  ecological  services  that  provide  these 
benefits  are  abundance  and  genetic  variance  (Bruggeman  et  al.  2005). 

The  goal  of  LEA  is  to  account  for  habitat  fragmentation  by  directing  the  trade  of  endangered 
species  habitat  toward  the  provision  of  the  ecological  services  provided  by  the  spatial  allocation 
of  habitat  in  which  the  organism  evolved.  Under  LEA,  the  financial  value  of  tradable  credits 
equals  the  in-kind  replacement  costs  of  ecological  services  that  are  differentially  affected  by 
changes  in  habitat  area  and  connectivity  (Bruggeman  et  al.  2005).  LEA  estimates  the 
equivalency  of  habitat  patches  traded  in  a  fragmented  landscape  based  on  changes  in  three 
ecological  services,  1)  abundance  and  genetic  variance  2)  within  and  3)  among  local  populations, 
measured  at  the  regional  scale.  Therefore,  a  patch  restored  is  equivalent  to  a  patch  lost  if  the 
alternative  landscape  structures  are  equivalent  in  all  three  ecological  services.  A  “LEA  credit” 
represents  the  marginal  contribution  a  change  in  landscape  structure  (i.e.,  taking  or  restoring  of 
breeding  or  dispersal  habitat)  makes  toward  moving  the  ecological  service  closer  to  or  away 
from  service  levels  observed  prior  to  habitat  loss  and  fragmentation. 

2.3  Method  for  Evaluating  Hypothetical  Trades  with  LEA  under  Alternative  Dispersal 
Models 

2.3.1  Step  1.  Constructing  Alternative  Landscapes 

In  order  to  estimate  LEA  credits,  ecological  service  flows  from  alternative  landscape 
configurations  must  be  projected.  LEA  credits  were  estimated  by  comparing  three  different 
landscape  configurations:  baseline,  withdrawal,  and  mitigation.  The  baseline  landscape  was 
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determined  using  MCBCL  RCW  habitat  restoration  program,  which  is  designed  to  meet  their 
contribution  to  the  Coastal  North  Carolina  Population  recovery  goal.  The  restoration  program  is 
focused  on  the  Main  Base  and  will  provide  a  highly  clumped  reference  landscape.  Levels  of 
genetic  variance  within  and  among  breeding  groups  expected  to  result  from  the  restoration 
program  were  used  as  the  baseline  levels  for  calculating  LEA  credits.  This  represents  an 
adjustment  from  the  “prior  to  habitat  loss  and  fragmentation”  condition  for  a  baseline  described 
in  Bruggeman  et  al.  (2005)  for  trading  among  private  landowners.  The  baseline  trajectory  of 
genetic  services  expected  by  the  restoration  plan  represents  the  improvement  in  rates  of 
recruitment,  migration,  and  local  extinction  that  result  from  meeting  yearly  recovery  objectives 
and  are  feasible  given  soil  types  and  existing  land  use.  Therefore,  habitat  trading  that  moves 
genetic  services  away  from  the  baseline  levels  has  resulted  in  moving  rates  of  migration,  drift, 
inbreeding,  and  local  extinction  away  from  levels  achievable  along  the  recovery  trajectory. 
Application  of  the  “prior  to  habitat  loss  and  fragmentation”  condition  applies  to  private 
landowners  because  they  trade  around  the  “take  standard”.  In  contrast,  Federal  landowners  trade 
around  a  “recovery  standard”.  Therefore,  baseline  service  trajectories  play  a  more  prominent 
role  in  credit  calculations  and  application  of  the  “prior  to  habitat  loss  and  fragmentation” 
condition  would  result  in  any  Federal  landowner  owing  credits  due  to  the  loss  of  endangered 
species  habitat  prior  to  the  ESA  being  promulgated. 

The  withdrawal  landscape  represents  the  landscape  condition  that  would  occur  due  to  training 
intensification  on  the  Main  Base,  which  would  result  in  taking  existing  clusters  or  not  restoring 
clusters  that  were  planned  to  contribute  to  recovery.  This  landscape  condition  can  be  compared 
with  a  landscape  in  which  new  clusters  are  established  outside  of  the  Main  Base  to  examine  the 
contribution  Encroachment  Partnering  makes  to  population  viability  in  light  of  expected  takes. 
This  landscape  will  be  referred  to  as  the  mitigation  landscape.  Below  we  describe  how  these 
landscapes  were  constructed  and  included  in  the  RCW  SEPMs. 

2. 3. 1.1  Baseline  Landscape.  MCBCL  plans  to  add  112  new  clusters  on  the  Main 
Base,  for  a  total  of  198  clusters  (Fig.  2.2),  to  meet  the  RCW  recovery  objective  of  173 
active  clusters.  Five  of  the  112  were  added  in  2006,  but  these  were  treated  as  planned 
clusters  because  our  project  began  before  monitoring  data  on  the  occupancy  of  these  new 
sites  were  available.  Using  MCBCL’s  “SppRePopPt”  GIS  layer,  we  determined  the 
cluster  locations  and  year  when  the  cluster  is  expected  to  be  made  available.  Within  the 
RCW  SEPMs,  we  added  the  cluster  at  the  time  and  location  planned  by  MCBCL.  We 
simulated  restoration  by  turning  immature  pine  and  mixed  pine-hardwood  pixels  into  the 
predominantly  pine  classification  (Table  1.1)  to  provide  100-hectares  of  foraging  habitat 
(USFWS  2003a).  We  later  learned  that  the  restoration  plan  present  in  the  GIS  layer  was 
a  rough  draft,  which  has  since  been  modified  to  restore  habitat  more  quickly  (Craig  Ten 
Brink,  MCBCL  ECON,  personal  communication). 

2.3. 1.2  Withdrawal  Landscape:  Proposed  Takes.  Due  to  training  intensification, 
MCBCL  has  identified  RCW  habitats  that  will  be  reduced  or  removed  in  order  to  achieve 
military  readiness.  MCBCL  provided  shapefiles  of  the  areas  to  be  affected  and  a  list  of 
territories  likely  to  be  lost  (Fig.  2.2).  Two  different  development  projects  are  examined 
(Table  2.1). 
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The  first  development  project  is  to  build  the  new  Marine  Special  Operations  Center 
(MARSOC),  in  the  southwest  corner  of  the  Main  Base,  beginning  in  2007  (Fig.  2.2). 
MARSOC  will  result  in  developing  88  hectares,  most  of  which  represents  immature 
forest,  previously  slated  for  restoration  of  Longleaf  pine  savannah  in  2020  to  contribute 
to  RCW  recovery  on  base.  Development  of  MARSOC  will  result  in  loss  of  one  planned 
cluster.  We  also  analyzed  a  worse  case  scenario  to  help  determine  how  sensitive  LEA  is 
for  analyzing  habitat  trades.  Therefore,  we  assumed  that  3  other  adjacent  clusters 
planned  for  installation  in  2020  would  not  be  included  in  the  simulation  (Table  2.1). 
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Figure  2.2.  Silhouette  of  Camp  Lejeune  used  to  highlight  proposed  changes  in  allocation  of  RCW  habitat.  Main  Base:  existing  RCW 
cluster  locations  are  indicated  with  black  #s  and  planned  clusters  in  white  #s;  development  of  MARSOC  and  the  G10  Area  shown 
with  stippled  areas.  Greater  Sandy  Run  Area:  12  potential  cluster  locations  are  numbered;  the  numbers  report  MCBCL’s  Forest 
Object  ID  for  stands  of  either  Loblolly  or  Longleaf  Pine  containing  at  least  2  pines  greater  than  10”  dbh.  Off-base:  parcels  purchased 
through  Encroachment  Partnering  (Beck,  Allen,  and  Lanier)  are  striped. 


34 


Final  Draft  v2 


Bruggeman,  Jones,  and  Wiegand 


The  second  development  project  is  the  installation  of  new  ranges  near  the  G10  area  on 
the  east  side  of  the  Main  Base  (Fig.  2.2).  Several  areas  of  existing  forest  stands  will  be 
developed,  totaling  232  hectares.  The  G10  development  is  expected  to  begin  in  2010  and 
impact  3  existing  clusters.  We  also  developed  a  worse  case  scenario  to  support  a 
sensitivity  analysis  for  LEA  by  assuming  that  5  more  clusters  near  these  development 
sites  would  be  lost  (Table  2.1). 

Table  2.1.  Description  of  proposed  takes  examined  with  Landscape  Equivalency  Analysis. 

_ Development  Project _ MARSOC _ G10  Area _ 

Hectares  affected  88  232 

Year  initiated  2007  2010 

Clusters  lost  -  Best  case  SB  HI  1  23,50,55 

Clusters  lost  -  Worst  case  SB  H8,  SB  H10,  SB  HI  1,  SB  H14  13,15,23,37,50,55,70,75 

Note:  planned  cluster  locations  are  identified  as  alpha-numeric,  while  existing  clusters  are  defined  numerically. 

2.3. 1.3  Mitigation  Landscapes:  Encroachment  Partnering  and  Restoration  of 
Greater  Sandy  Ran  Area.  Three  different  mitigation  landscapes  were  evaluated  with 
LEA.  We  evaluated  the  conservation  benefit  of  restoring  three  of  the  five  parcels  already 
purchased  through  MCBCL’s  Encroachment  Partnering  program.  These  parcels  were 
chosen  due  to  their  close  adjacency  to  the  Main  Base  and  potential  for  allowing  habitat 
trades  without  exacerbating  fragmentation.  The  parcels  included  the  Beck,  Allen,  and 
Lanier  tracts  (Fig  2.2).  We  also  considered  the  conservation  benefit  of  restoring  the 
Greater  Sandy  Run  Area  (GSRA),  which  is  adjacent  to  the  Main  Base  and  solely  owned 
by  MCBCL. 

Mitigation  Choice  1 :  Restore  only  the  Beck  and  Allen  tracts.  The  Beck  and  Allen  tracts 
are  adjacent  to  each  other  and  lie  south  of  MCBCL  and  west  of  the  New  River  (Fig.  2.2). 
The  Beck  tract  is  976  hectares  and  was  acquired  in  cooperation  with  NC  Wildlife 
Resources  Commission  (NCWRC).  NCWRC  contributed  $3,396,000  and  MCBCL 
contributed  $5,236,000,  for  a  total  cost  of  $8,632,000.  The  Allen  tract,  sometimes 
referred  to  as  Testamentary,  is  68  hectares,  and  was  acquired  with  cooperation  with  North 
Carolina  Coastal  Land  Trust  (NCCLT).  MCBCL  and  NCCLT  both  contributed 
$1,732,500  (total  cost  =  $3,465,000).  In  total,  the  Beck  and  Allen  tracts  comprise  1,044- 
ha. 

The  vegetative  composition  of  these  lands  has  not  yet  been  quantified.  Based  on  an  aerial 
photograph,  these  tracts  contain  mixed  pine -hardwood  forests  and  bottomland  hardwood 
forests.  MCBCL  estimates  that  8  RCW  clusters  can  be  located  on  the  Beck  tract  and  one 
on  the  Allen  tract  in  about  30  years,  with  an  average  territory  size  of  100-ha  (USFWS 
2003a).  A  GIS  shape  file  with  property  boundaries  was  provided.  Based  on  soil 
classification,  MCBCL  provided  a  photo  of  the  approximate  areas  where  Longleaf  could 
be  restored.  We  used  this  information  to  estimate  cluster  centers  and  integrated  the 
current  and  expected  habitat  conditions,  derived  from  viewing  an  aerial  photograph  and 
soils  map,  into  the  landscape  simulation  model. 
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Mitigation  Choice  2:  Restore  only  the  Lanier  tract.  The  Lanier  tract  is  a  177-hectare 
parcel  and  lies  in  between  the  Main  Base  and  GSRA  to  the  north  (Fig.  2.2).  It  was 
purchased  with  a  cooperative  agreement  between  MCBCL  and  NCWRC  for  a  total  cost 
of  $2,032,500.  Based  on  aerial  photograph,  the  Lanier  tract  is  comprised  mostly  of 
mixed  pine-hardwood  forests  and,  secondarily,  bottomland  hardwood  forests.  We 
estimated  that  two  clusters  could  be  located  on  the  Lanier  tract  in  30  years. 

Mitigation  Choice  3:  Restore  Beck,  Allen,  Lanier,  and  Greater  Sandy  Run  Area. 

We  did  not  consider  restoration  of  GSRA  without  the  other  three  parcels  because  the 
Beck,  Allen,  and  Lanier  tracts  will  increase  the  connectivity  of  clusters  in  GSRA  with 
those  on  the  Main  Base.  MCBCL  has  initiated  the  restoration  of  Longleaf  Pine  in  GSRA 
with  the  hopes  of  better  connecting  their  RCW  population  with  birds  in  the  Holly  Shelter 
Games  Lands  and  other  parcels  under  consideration  for  the  Encroachment  Partnering 
program  (Fig  2. 1).  GSRA  was  acquired  in  1991  by  MCBCL  from  International  Paper 
primarily  for  military  training. 

Placement  of  new  clusters  in  the  Greater  Sandy  Run  Area  was  based  on  the  best  available 
data  from  MCBCL’s  Forestry  Coverage  GIS  Layer  and  the  RCW  Foraging  Habitat 
Guidelines  (USFWS  2003a).  The  Foraging  Habitat  Guidelines  require  information  on  the 
basal  area  of  specific  size  class  distributions  of  pine.  However,  these  data  were  often  not 
available  for  the  older  stands  of  Longleaf  and  Loblolly  pines  in  the  GSRA.  One  basic 
requirement  of  the  guidelines  is  that  good  quality  foraging  habitat  should  have  40  pines 
per  hectare  (or  16  pines  per  acre)  of  pines  at  least  36  cm  dbh  and  cover  5 1  hectares. 
Available  data  from  Forestry  Coverage  GIS  Layer  allowing  us  to  infer  these 
characteristics  included  birth  year  of  a  representative  pine,  total  pine  basal  area,  primary 
species,  number  of  pines  greater  than  25  cm  dbh,  and  size  of  stand  (Table  2.2). 

We  identified  12  forest  stands  in  GSRA  dominated  by  Longleaf  or  Loblolly  Pine  and 
containing  at  least  2  pines  greater  than  25  cm  dbh  (Table  2.2;  Fig  2.2).  These  stands 
tended  to  be  small  in  size  and  surrounded  by  Longleaf  and  Loblolly  stands  that  were 
much  younger.  Therefore,  we  assumed  that  stands  with  a  birth  year  less  than  1960, 
would  not  be  suitable  for  both  reproductive  and  foraging  habitat  for  another  30  years. 

For  the  remaining  six  stands,  we  assumed  that  they  would  be  suitable  within  50  years. 
Restoration  was  simulated  by  turning  immature  pine  and  mixed  pine-hardwood  pixels 
within  800  m  of  cluster  center  into  the  predominantly  pine  classification  (Table  1.1)  to 
provide  100-hectares  of  foraging  habitat  (USFWS  2003a). 
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Table  2.2.  Forest  attributes  describing  stands  in  Greater  Sandy  Run  Area  that  may  be  suitable 
locations  for  RCW  nesting  cavities. 

Number  of 


Year  cluster 
include  in 
Simulation 

Stand 

Object 

ID 

Birth 
year  of 
Stand 

Pine 

basal 

Area 

Hardwood 

basal 

Area 

pines 

greater  than 
25  cm  dbh 

Primary 

Species 

Hectares 

320 

1939 

83 

0 

2 

Longleaf 

2 

2035 

2306 

1948 

70 

0 

4 

Loblolly 

4 

2035 

2691 

1948 

76 

0 

2 

Loblolly 

20 

2035 

3327 

1948 

70 

0 

2 

Longleaf 

4 

2035 

501 

1953 

68 

0 

3 

Longleaf 

4.5 

2035 

1070 

1956 

95 

0 

3 

Longleaf 

3.6 

2045 

69 

1968 

75 

0 

2 

Loblolly 

18 

2045 

261 

1968 

99 

0 

4 

Loblolly 

28 

2045 

215 

1970 

100 

0 

6 

Loblolly 

15 

2055 

1819 

1973 

80 

0 

4 

Loblolly 

1.2 

2055 

785 

1976 

150 

0 

2 

Loblolly 

53 

2836 

1980 

90 

0 

2 

Loblolly 

22 

2.3.2  Step  2.  Projecting  Ecological  Services  under  Uncertainty 

In  order  to  predict  ecological  services  in  the  future  for  each  landscape  configuration,  we  used  the 
Letcher  model  and  eight  separate  models  of  the  new  HCRW  SEPM  identified  in  Phase  1 .  Each 
of  the  eight  HCRW  models  represents  a  different  parameterization  that  passed  through  all  five 
pattern  fulfillment  criteria  using  Pattern-Oriented  Modeling.  Using  multiple  models  allows  us  to 
incorporate  how  uncertainty  regarding  a  species’  behavior  and  demography  might  affect  decision 
making. 

Each  SEPM  was  initiated  with  the  population  structure  observed  on  MCBCL  in  2005,  including 
the  number  of  observed  breeders,  helpers,  fledglings,  and  floaters.  Initial  population  size  was 
316  individuals.  All  models  included  an  infinite  alleles  model  of  genetic  variance  (Lacy  2000). 
The  infinite  alleles  model  assumes  that  every  founding  individual  is  heterozygous  and  contains 
two  unique  alleles  (i.e.,  total  alleles  =  2  x  number  of  breeders  in  the  founding  population).  We 
compiled  a  pedigree  using  MCBCL’s  RCW  monitoring  data  provided  by  Jeffrey  Walters.  Using 
the  pedigree  we  simulated  the  transmission  of  alleles  from  parents  to  offspring  assuming 
Mendelian  inheritance  (i.e.,  their  offspring  had  an  equal  probability  of  inheriting  each  of  the  two 
alleles)  for  the  years  1986-2005,  at  the  start  of  each  simulation.  Simulating  Mendelian 
inheritance  under  an  infinite  alleles  model  allows  us  to  estimate  the  spatial  distribution  of  genetic 
variance  in  the  initial  population. 

Population  dynamics  were  projected  100  years  into  the  future.  At  the  start  of  each  year,  eight 
immigrants,  or  floaters,  were  added.  The  average  number  of  immigrants  observed  each  year  in 
MCBCL  over  the  past  five  years  was  8.2.  We  assumed  that  each  immigrant  was  unrelated  to 
each  other  and  all  individuals  on  MCBCL.  Therefore,  each  immigrant  was  heterozygous, 
containing  two  unique  alleles.  Immigrants  were  added  within  the  simulation  to  parallel  the 
manner  in  which  Jeffrey  Walters’  field  crews  first  observe  new  floaters.  Floaters  of  unknown 
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parentage  are  observed  in  the  spring  near  a  cluster.  Therefore,  we  placed  each  new  floater  at  a 
location  randomly  selected  from  all  possible  clusters  available  in  a  year,  whether  occupied  or 
not.  At  that  point  the  immigrant  can  either  compete  for  a  breeding  vacancy  within  their 
perceptual  distance  or  begin  their  stepwise  search  for  breeding  vacancies  elsewhere. 

In  order  to  estimate  LEA  credits,  ecological  services  in  the  baseline,  withdrawal,  and  mitigation 
landscapes  were  calculated.  We  estimated  abundance  as  the  number  of  PBGs  due  to  the 
biological  and  regulatory  significance  of  this  estimate  (USFWS  2003a).  Recall  that  MCBCL’s 
contribution  to  recovery  will  be  estimated  based  on  number  of  “active  clusters”.  Due  to  the 
difficulty  of  estimating  when  neighboring  groups  merely  capture  a  cluster  but  do  not  use  it  for 
breeding  within  a  computer  simulation,  we  used  PBGs  instead.  The  expected  apportionment  of 
genetic  variance  within  and  among  breeding  groups  under  the  assumption  of  Hardy- Weinberg 
Equilibrium  was  estimated  by  the  method  of  Nei  (1973;  as  described  in  Appendix  A).  Nei’s 
theory  relates  the  total  expected  heterozygosity  in  the  population  (H  r)  to  the  average  expected 

heterozygosity  within  breeding  groups  (H§)  and  average  genetic  divergence  among  breeding 
groups  (Dst),  Ht  =  Hs  +  Dst-  The  SEPMs  include  stochastic  effects  when  detennining  the 
number  of  fledglings  produced,  individual  survival,  and  dispersal  behaviors  as  previously 
discussed  in  Phase  1.  Therefore,  we  took  an  average  for  PBG,  Hs,  and  DSt  at  each  year  from  100 
replicate  simulations. 

2.3.3  Step  3.  Estimate  Spatially-explicit  Debits  and  Credits  with  LEA 

LEA  estimates  debits  and  credits  using  a  time-integrated  estimate  of  proportional  differences  in 
ecological  services  expected  under  alternative  landscape  configurations,  called  Landscape 
Service  Years  (LSYs).  LSYs  provide  an  estimate  of  the  change  in  landscape  function  due  to 
changing  habitat  area  and/or  connectivity. 

The  debit  in  abundance  is  estimated  as  the  proportional  change  in  PBG  due  to  habitat  removal 
relative  to  number  of  PBG  expected  under  the  recovery  trajectory.  The  debit  in  abundance 
services  due  to  habitat  loss  (LSYdPBG)  represents  the  marginal  conservation  value  of  the  clusters 
lost,  given  their  spatial  context  within  the  population. 


Where  the  expected  loss  of  services  due  to  habitat  loss  begins  at  year  W,  rtPBG  is  the  expected 
number  of  PBGs  at  year  t  under  the  USFWS  recovery  trajectory,  jtPBG  is  the  expected  number  of 
PBGs  at  year  t  given  the  loss  of  habitat  (Table  2.1). 


Next  we  estimate  the  increase  in  services  due  to  restoring  habitat  elsewhere.  The  expected 
number  of  abundance  credits  generated  by  the  mitigation  landscape  (LSYCPBG)  is  estimated  as: 


(eqn  2.2) 


,  PBG 

where  mt 


is  the  expected  number  of  PBGs  at  year  t  given  the  loss  of  habitat  on  the  Main  Base 


and  restoration  of  habitat  outside  of  the  Main  Base  (i.e.,  service  trajectory  from  mitigation 
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landscape).  If  (LSYcPBCt  -  LSYdPBG)  is  greater  than  or  equal  to  zero,  then  off-site  restoration 
provides  enough  abundance  credits  to  justify  the  habitat  trade. 

Calculating  tradable  credits  associated  with  changes  in  genetic  variance  is  more  complex.  The 
management  goal  is  to  approximate  ecological  services  provided  by  the  distribution  of  habitat  in 
which  the  organism  evolved  (Meffe  1996;  Moritz  2002).  Greater  genetic  diversity  within  a 
breeding  group  or  greater  genetic  divergence  among  breeding  groups  is  not  always  better  for 
sustainability  (Bouzat  2001).  The  levels  of  genetic  service  provided  by  the  baseline  landscape 
will  be  used  to  direct  trading  toward  this  goal.  The  spatial  distribution  of  the  clusters  in  the 
baseline  landscape  have  been  planned  to  be  as  clumped  in  distribution  as  vegetative  composition 
will  allow. 


As  estimates  of  genetic  variance  within  and  among  breeding  groups  move  farther  away  from 
baseline  levels  due  to  losing  a  patch  or  connectivity,  the  larger  the  debit  in  genetic  services  that 
accrues.  To  account  for  changes  in  genetic  diversity  either  greater  than  or  less  than  baseline 
levels,  LSYd  uses  the  absolute  distance  between  service  trajectories: 


w 


100 

=  z 

t=W 


b^- 


Jt 


V 


(eqn  2.3) 


where  G  is  the  genetic  variance  component  estimated  (Hs  or  DST),  bt  1  is  the  expected  level  of 

p 

genetic  variance  at  year  t  provided  by  the  baseline  landscape,  and  jt  is  the  expected  level  of 
genetic  variance  at  year  t  provided  by  the  withdrawal  landscape. 


To  estimate  the  contribution  restoring  habitat  off-base  would  make  toward  maintaining  the 

p 

spatial  allocation  of  genetic  variance  observed  in  the  baseline  landscape,  LSYC  is  estimated  as: 


LSYr 


100 
=  Z 
J  t=W 


bt^~ 


Jt 


100 

-  z 

t=W 


b y- 


mt 


(eqn  2.4) 


p 

where  mt  J  is  the  expected  level  of  genetic  variance  at  year  t  provided  by  the  mitigation 
landscape,  which  includes  the  loss  of  clusters  on  the  Main  Base.  LSYC  is  large  when  mitigation 
moves  genetic  services  close  to  baseline  (i.e.,  value  on  right  hand  side  of  minus  sign  equals 
zero).  The  best  outcome  is  that  LSYC  J  equals  LSYd  ,  meaning  that  clusters  on-base  can  be 
traded  with  clusters  outside  of  the  Main  Base  without  exacerbating  the  effects  of  habitat 
fragmentation.  In  other  words,  restoration  has  been  provided  so  that  the  baseline  and  mitigation 
landscapes  provide  equivalent  genetic  services.  IfLSYc  is  positive  but  less  than  LSYd  ,  then 
mitigation  off-base  did  help  ameliorate  the  effects  of  habitat  fragmentation  caused  by  the  take, 
but  rates  of  recruitment  and  migration  in  the  baseline  and  mitigation  landscapes  are  not 
equivalent.  If  LSYC 1  equals  zero  then  the  withdrawal  and  mitigation  level  genetic  services  are 
equally  distant  from  baseline  levels  and  mitigation  did  not  help  to  reduce  the  effects  of 
fragmentation.  If  LSYC  is  negative  then  mitigation  moves  genetic  services  farther  away  from 
baseline  than  observed  under  the  withdrawal  landscape  indicating  that  off-base  mitigation  has 
exacerbated  functional  fragmentation. 
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2.4  Results  and  Discussion 

2.4.1  Mitigation  Choice  1:  Restore  Only  Beck  and  Allen  Tracts 

First,  we  discuss  the  extent  to  which  restoring  9  clusters  on  the  Beck  and  Allen  tracts  will 
mitigate  the  loss  of  clusters  on  the  Main  Base  (Table  2. 1),  including  the  secondary  effects  of 
changes  in  habitat  connectivity.  The  temporal  trajectories  of  ecological  services  are  illustrated  in 
Figures  2.3  to  2.5.  Recall  that  MCBCL  is  required  to  reach  its  recovery  goal  by  2025.  However, 
habitat  loss  begins  in  2007  while  habitat  restoration  outside  of  the  Main  Base  does  not  begin 
until  2035. 

We  observed  that  the  withdrawal  landscape  will  still  provide  enough  clusters  to  exceed  the 
recovery  goal  by  2025  (Fig  2.3),  thanks  to  the  aggressive  RCW  habitat  restoration  plan  on  Main 
Base.  However,  the  rate  of  recovery  in  the  withdrawal  landscape  is  a  little  slower  when 
compared  to  the  baseline  landscape,  especially  in  the  Letcher  model.  A  high  rate  of  habitat 
occupancy  was  observed  for  the  Letcher  model  and  all  8  parameterizations  of  the  HCRW  SEPM. 
Therefore,  the  debit  in  services  estimated  by  comparing  the  PBG-trajectory  under  the  withdrawal 
landscape  to  the  PBG-trajectory  given  the  recovery  goal  (eqn  2.1)  resulted  in  a  negative  LSYdPBG 
value,  because  the  j -trajectory  exceeds,  on  average,  the  r-trajectory  from  the  time  the  loss  of 
habitat  begins  (W=2007)  (Table  2.3).  In  other  words,  despite  the  loss  of  habitat  area  in  both  the 
best  and  worst  case  scenarios,  all  models  predict  that  more  than  enough  clusters  are  provided  and 
occupied  by  RCWs  on  the  Main  Base  to  meet  MCBCL’s  contribution  to  recovery  by  2025. 

The  conservation  benefit  of  establishing  new  clusters  on  the  Beck  and  Allen  tracts  in  2035, 
estimated  by  the  m-trajectory,  indicates  occupancy  of  these  new  clusters  shortly  after  they  were 
made  available  in  2035  (Fig  2.3).  Therefore,  LSYCPBG  values  indicate  a  further  accrual  of  PBG 
credits  relative  to  the  recovery  goal,  for  all  models  (eqn.  2.2;  Fig  2.3). 

The  differences  in  LEA  credits  between  the  new  RCW  HCRW  and  the  Letcher  SEPMs  are  not 
large  for  PBGs  (Table  2.3).  Both  indicate  that  more  credits  remain  under  the  best  case 
development  scenario.  The  Letcher  model  predicts  that  habitat  loss  will  have  a  slightly  worse 
effect  on  PBGs  than  the  HCRW  models.  Further,  the  Letcher  model  predicts  that  habitat 
restoration  will  have  a  slightly  stronger  benefit  for  habitat  occupancy  than  predicted  by  HCRW 
models.  The  net  result  is  that  both  the  Letcher  and  HCRW  models  predict  roughly  the  same 
number  of  LSYPBG  credits  will  remain  after  the  trade. 

Proportional  changes  in  genetic  services  due  to  the  loss  of  clusters  on  the  Main  Base  were  small 
(Table  2.3;  Figs.  2.4  and  2.5).  Comparing  the  average  expected  heterozygosity  within  breeding 
groups  (Hs)  in  the  baseline  landscape  to  that  observed  in  the  withdrawal  landscape  (i.e.,  eqn.  2.3 
reporting  LSYdHs)  suggested  that  the  loss  of  clusters  on  the  Main  Base  may  only  slightly  change 
rates  of  genetic  drift  and  inbreeding  (Table  2.3;  Fig  2.4).  The  Letcher  model  predicts  LSYdHs- 
values  that  fall  within  the  range  expected  from  the  eight  HCRW  models.  The  rather  small  debit 
for  Hs  is  not  surprising  given  the  density  of  RCW  clusters  that  remain  on  the  Main  Base  after  the 
MARSOC  and  G10  development  projects,  and  the  number  of  immigrants  expected  to  arrive  each 
year  to  MCBCL. 
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When  including  the  comparison  of  the  baseline  landscape  to  the  mitigation  landscape  (eqn.  2.4), 
differences  between  the  Letcher  and  HCRW  models  emerge.  The  Letcher  model  predicts  that 
off-base  restoration  will  help  move  Hs  closer  to  baseline  levels  (i.e.,  LSYC  values  are  positive). 
However,  LSYcHs  values  are  not  as  large  as  LSYdHs  values,  so  a  small  habitat  fragmentation 
effect  may  remain  despite  mitigation  efforts.  In  contrast,  the  HCRW  models  predict  that  off-base 
mitigation  will  move  Hs  farther  away  from  baseline  (i.e.,  LSYcHs  values  are  negative). 

Therefore,  off-base  mitigation  has  increased  the  effects  of  habitat  fragmentation  compared  to  the 
withdrawal  landscape.  The  negative  values  for  LSYcHs  are  small,  but  they  suggest  that  the  9  off- 
base  clusters  do  not  make  an  equivalent  contribution  to  rates  of  recruitment  as  those  that  may  be 
lost  on  the  Main  Base.  Recall  the  Letcher  model  does  not  allow  land  cover  to  affect  dispersal 
behaviors.  Therefore,  breeding  vacancies  may  be  filled  at  a  greater  rate  in  the  Letcher  model. 

For  average  genetic  divergence  among  breeding  groups  (Dst),  we  observed  a  small  debit  due  to 
the  loss  of  clusters  on  the  Main  Base  (LSYd°st;  Fig.  2.5).  The  loss  of  these  clusters  often  resulted 
in  an  increase  in  genetic  differences  among  groups.  In  other  words,  the  loss  of  the  clusters  on  the 
Main  Base  likely  changed  patterns  of  migration  among  groups,  leading  to  an  increase  in  the 
degree  of  genetic  differentiation  among  groups.  When  comparing  DST  values  observed  in  the 
baseline,  withdrawal,  and  mitigation  landscapes  (eqn.  2.4),  a  small  increase  in  habitat 
fragmentation  effects  is  observed  in  all  models  (Table  2.3).  Therefore,  the  9  clusters  restored 
off-base  do  not  make  equivalent  contributions  to  migration  among  clusters  as  those  that  may  be 
lost  on  the  Main  Base.  However,  the  net  increase  in  habitat  fragmentation  effects  is  small  (Fig 
2.5).  The  Letcher  model  predicts  the  smallest  increase  in  genetic  divergence  because  its 
dispersal  rules  allow  for  relatively  higher  rates  of  gene  flow  among  distant  clusters. 
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Table  2.3.  Landscape  Service  Years  (LSYs)  that  evaluate  the  ability  of  restoring  9  clusters  on  the  Beck  and  Allen  tracts  to  offset 
habitat  loss  and  fragmentation  effects  resulting  from  the  loss  of  4  clusters  (best  case)  or  12  clusters  (worst  case)  on  the  Main  Base. 

Landscape  Service  Years  (LSYs) 

PBG  Credits  4Hs  Credits  4DSt  Credits 


2LSYdPBG 

LSYCPBG 

remaining 

LSYdHs 

3LSYcHs 

remaining 

LSYdDst 

3LSYcDst 

remaining 

XHCRW  SEPM 

Best  Case 

minimum 

-8.38 

3.17 

9.83 

0.084 

-0.212 

-0.297 

0.203 

-0.498 

-0.705 

median 

-7.60 

3.34 

11.00 

0.097 

-0.146 

-0.251 

0.228 

-0.392 

-0.635 

maximum 

-6.54 

3.55 

11.92 

0.154 

-0.080 

-0.195 

0.288 

-0.259 

-0.514 

Letcher  SEPM 

-6.63 

3.60 

10.23 

0.096 

0.004 

-0.092 

0.173 

-0.119 

-0.293 

!HCRW  SEPM 

Worst  Case 

minimum 

-4.13 

2.83 

5.26 

0.068 

-0.271 

-0.375 

0.183 

-0.581 

-0.900 

median 

-3.31 

3.21 

6.52 

0.104 

-0.141 

-0.236 

0.251 

-0.343 

-0.599 

maximum 

-2.43 

3.38 

7.51 

0.126 

-0.068 

-0.161 

0.327 

-0.195 

-0.429 

Letcher  SEPM _ -2.46  3.43 _ 5.89  0.096  0.022  -0.074  0.226  -0.035 _ -0.262 

1  Range  of  values  taken  from  results  of  running  eight  different  parameterizations  of  the  RCW  HCRW  SEPM  derived  from  Phase  1 . 

2  Proportional  loss  of  ecological  services  relative  to  the  recovery  goal  (eqn  2.1).  If  negative,  a  gain  in  ecosystem  services  has 
occurred  above  that  required  by  the  recovery  criteria. 

3  Mitigation  exacerbated  the  effects  of  fragmentation  when  values  are  negative  (eqn  2.4). 

4  When  values  are  close  to  zero,  habitat  trading  did  not  exacerbate  the  effects  of  fragmentation. 


42 


Final  Draft  v2 


Bruggeman,  Jones,  and  Wiegand 


CL 
LU 
C 0 


0 

-£= 

o 

0 


00 

a" 

w  o 
^  fo 

O  N 
*'§ 
^  0 
X  £ 
O  2 

X  CD 
CL 


Best  Case  Scenario  Worst  Case  Scenario 


200 

. AAAAaaAAAAAAAAAAAAAAAAAAAAAAAAAAA 

200 

•  180- 

•  180- 

=  160 

JT 

=  160 

o 

Jr 

o 

5  140 

5  140 

/f 

c  120 

/A 

c  120 

/A 

TJ 

/JF 

■o 

S  ioo 

b  -  baseline 

g  100 

jfccfc)  h  -  h3s;filinfi 

i- 

B 

co 

o 

Sf  O  j  -  withdrawal 

B 

00 

o 

Jp  o  j  -  withdrawal 

~  60 

t  Am-  mitigation 

=  60 

f  Am-  mitigation 

ote 

O 

^  r  -  Recovery  Goal 

ote 

4^ 

O 

^  r  -  Recovery  Goal 

Q- 

Q. 

20 

■  Observed 

20 

■  Observed 

1111111 

1  1  1  1  1  1  1 

1995  2010  2025  2040  2055  2070  2085  2100 

1995  2010  2025  2040  2055  2070  2085  2100 

200 

^....AAAAAAAAAAAAAAAAAAAAAAAAAA 

200 

•  I80- 

•  180- 

=  160 

<T 

=  160 

o 

A 

o 

r 

5  140 

/I 

5  140 

/I 

.WJ 

jTJ 

c  120 

c  120 

T3 

TJ 

$  100 

b  -  baseline 

g  100 

h  -  basfilinft 

1- 

1- 

B 

00 

o 

O  j  -  withdrawal 

B 

00 

o 

J?  o  j  -  withdrawal 

03 

~  60 

f  Am-  mitigation 

~  60 

f  Am-  mitigation 

f  40 

^  *  r  -  Recovery  Goal 

%  40 

^  r  -  Recoery  Goal 

“■  20 

■  Observed 

“■  20 

■  Observed 

1111111 

1111111 

1995  2010  2025  2040  2055  2070  2085  2100 

1995  2010  2025  2040  2055  2070  2085  2100 

Figure  2.3.  Trajectories  of  Potential  Breeding  Groups  used  to  estimate  LEA  credits  derived  from  the  original  Letcher  SEPM  and  one 
parameterization  of  the  Habitat-correlated  Random  Walk  SEPM.  Baseline  landscape  includes  the  addition  of  1 12  new  clusters  by 
2025.  Withdrawal  landscape  represents  the  baseline  landscape  minus  the  loss  of  clusters  due  to  MARSOC  (t=2007)  and  G10  (t=20 10) 
development.  The  best  and  worst  case  scenarios  are  defined  in  Table  2.1.  The  mitigation  landscape  represents  the  withdrawal 
landscape  plus  the  addition  of  clusters  on  the  Beck  and  Allen  tracts  (t=2035)  thanks  to  Encroachment  Partnering.  The  recovery 
trajectory  represents  USFWS’s  expectation  for  MCBCL’s  contribution  to  recovery  of  the  Coastal  North  Carolina  RCW  Population. 
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Best  Case  Scenario 


Worst  Case  Scenario 
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Figure  2.4.  Trajectories  of  average  expected  heterozygosity  within  breeding  groups  (Hs)  used  to  estimate  LEA  credits  derived  from 
the  original  Letcher  SEPM  and  one  parameterization  of  the  Habitat-correlated  Random  Walk  SEPM.  Parameterization  77481  is 
displayed  because  the  biggest  change  in  genetic  diversity  was  observed  in  this  simulation.  Baseline  landscape  includes  the  addition  of 
112  new  clusters  by  2025.  Withdrawal  landscape  represents  the  baseline  landscape  minus  the  loss  of  clusters  due  to  MARSOC 
(t=2007)  and  G10  (t=2010)  development.  The  best  and  worst  case  scenarios  are  defined  in  Table  2.1.  The  mitigation  landscape 
represents  the  withdrawal  landscape  plus  the  addition  of  clusters  on  the  Beck  and  Allen  tracts  (t=2035)  thanks  to  Encroachment 
Partnering. 
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Best  Case  Scenario 


Worst  Case  Scenario 


Figure  2.5.  Trajectories  of  average  genetic  divergence  among  breeding  groups  (Dst)  used  to  estimate  LEA  credits  derived  from  the 
original  Letcher  SEPM  and  one  parameterization  of  the  Habitat-correlated  Random  Walk  SEPM.  Parameterization  77481  is  displayed 
because  the  biggest  change  in  genetic  diversity  was  observed  in  this  simulation.  Baseline  landscape  includes  the  addition  of  1 12  new 
clusters  by  2025.  Withdrawal  landscape  represents  the  baseline  landscape  minus  the  loss  of  clusters  due  to  MARSOC  (t=2007)  and 
G10  (t=2010)  development.  The  best  and  worst  case  scenarios  are  defined  in  Table  2.1.  The  mitigation  landscape  represents  the 
withdrawal  landscape  plus  the  addition  of  clusters  on  the  Beck  and  Allen  tracts  (t=2035)  thanks  to  Encroachment  Partnering. 
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2,4.2  Mitigation  Choice  2:  Restore  only  Lanier 

Next  we  contrast  the  conservation  benefit  of  restoring  only  the  Beck  and  Allen  tracts  with  the 
choice  to  restore  only  the  Lanier  tract.  The  trades  are  compared  based  on  the  number  of  LSYs 
remaining  after  a  trade  (i.e.,  LSYd  -  LSYC),  which  represents  the  net  gain  or  loss  of  ecological 
services  due  to  a  change  in  landscape  pattern. 

As  expected,  restoring  only  2  clusters  on  the  Lanier  tract  will  leave  fewer  LEA  credits  for  PBGs 
than  would  be  possible  if  the  Beck  and  Allen  tracts  were  restored  (Fig  2.6).  These  results  do 
indicate  that  the  Lanier  clusters  also  have  a  high  probability  of  maintaining  breeding  groups 
(LSYcpbg  >  0  for  all  models).  For  all  off-base  habitat  trading,  the  number  of  FSYs  remaining  for 
PBGs  reported  by  the  Fetcher  model  is  within  the  range  of  values  expected  from  the  eight 
HCRW  models  (Fig  2.6). 

Changes  in  genetic  services  indicate  that  habitat  fragmentation  effects  are  smaller  when  restoring 
the  Fanier  tract  than  the  Beck  and  Allen  tracts  (Fig.  2.6).  In  other  words,  restoring  only  2 
clusters  on  Fanier  provided  more  equivalent  contributions  to  the  rates  of  recruitment  and 
migration,  compared  to  those  lost  on  Main  Base,  than  the  9  clusters  that  could  be  restored  on  the 
Beck  and  Allen  tracts.  Therefore,  while  the  Lanier  tract  represents  a  smaller  increase  in  habitat 
area,  it  provides  clusters  that  are  more  strongly  connected  with  other  clusters  in  the  landscape. 
The  Letcher  model  suggests  that  the  apportionment  of  genetic  variance  is  fairly  insensitive  to 
changes  in  landscape  structure  that  result  from  trading  (Fig  2.6).  In  contrast,  the  HCRW  models 
suggest  that  restoration  of  Lanier  will  result  in  smaller  changes  in  population  genetic  structure, 
compared  to  the  restoration  of  the  Beck  and  Allen  tracts.  The  changes  in  genetic  services  due  to 
mitigation  is  often  very  small  (i.e.,  the  largest  possible  change  in  genetic  structure  for  these 
comparisons  is  displayed  in  Figs.  2.4  and  2.5). 


■  HCRW  min 


j  =  best  case,  j  =  worst  case,  j  =  best  case,  j  =  worst  case, 
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j  =  best  case,  j  =  worst  case,  j  =  best  case,  j  =  worst  case, 
m  =  Beck-Alien  m  =  Beck-Alien  m  =  Lanier  m  =  Lanier 
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Figure  2.6.  LEA  credits  remaining  after 
trading  on-base  clusters  for  those  on  the 
Beck  and  Allen  tracts  or  for  those  on  the 
Lanier  tract.  Results  for  the  minimum 
and  maximum  number  of  LSYs  from  the 
eight  HCRW  SEPMs  and  the  Letcher 
SEPM  are  used  to  evaluate  changes  in 
RCW  population  dynamics. 
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2.4.3  Mitigation  Choice  3:  Restore  Beck,  Allen,  Lanier  and  GSRA 

We  examined  the  conservation  benefit  of  restoring  Beck,  Allen,  Lanier,  and  GSRA  for  the  best 
case  and  worst  case  development  scenarios.  Restoring  all  possible  areas  outside  of  the  Main 
Base  should  help  increase  connectivity  of  RCWs  across  the  Onslow  Bight  Landscape.  However, 
loss  of  forested  areas  due  to  economic  growth  is  occurring  rapidly  across  the  Onslow  Bight 
Landscape.  Therefore,  we  included  a  new  question:  how  might  further  development  in  the 
Onslow  Bight  Landscape  affect  future  trading  options  for  MCBCL?  We  tested  this  by  varying 
the  number  of  immigrants  added  per  year  to  the  simulation.  Recall  that  eight  floaters  of 
unknown  origin  are  added  to  the  landscape  at  eight  randomly  selected  clusters  each  year  to 
parallel  field  observations.  To  simulate  how  development  may  affect  landscape  function  we 
evaluated  trades  assuming  8,  2,  and  0  immigrants  per  year. 

First,  we  compare  the  results  of  the  three  mitigation  choices  by  assuming  eight  immigrants  per 
year.  The  restoration  of  Beck,  Allen,  Lanier,  and  GSRA  provides  an  increase  in  credits 
remaining  for  PBGs  (Fig.  2.7),  compared  to  the  decisions  to  develop  only  the  Beck  and  Allen 
tracts  or  only  the  Lanier  tract  (Fig.  2.6).  However,  comparing  Figures  2.6  and  2.7,  we  see  that 
the  increase  in  credits  is  not  linear.  In  other  words,  if  the  Beck,  Allen,  and  Lanier  tracts  achieved 
the  same  level  of  habitat  occupancy  when  all  mitigation  sites  were  restored,  then  we  would 
expect  that  we  could  simply  stack  the  histograms  from  Figures  2.6  and  2.7.  In  the  best  case 
scenario,  if  we  stacked  the  histograms,  we  would  get  18-22  LSYs  remaining  for  PBGs  instead  of 
12  to  16  (Fig  2.7).  Therefore,  providing  more  clusters  that  are  dispersed  may  help  decrease 
occupancy  previously  observed  under  the  Beck  and  Allen  or  Lanier  only  mitigation  scenarios 
(Pulliam  1988).  When  immigration  into  MCBCL  is  decreased  (Fig.  2.7),  a  small  decline  in  the 
number  of  LSYs  remaining  for  PBGs  results  in  both  the  best  and  worst  case  scenarios. 

By  restoring  all  sites  outside  of  the  Main  Base,  we  observed  that  habitat  fragmentation  effects 
increased  notably  in  the  HCRW  models.  Assuming  eight  immigrants  per  year,  the  decrease  in 
credits  remaining  observed  for  Hs  under  mitigation  choice  3  (Fig  2.7),  compared  to  the  other 
mitigation  choices  (Fig.  2.6),  indicate  that  the  new,  dispersed  clusters  do  not  make  equivalent 
contributions  to  rates  of  recruitment  as  those  that  may  be  lost  from  the  Main  Base.  The  Letcher 
model  also  indicates  increased  fragmentation  effects  due  to  mitigation  choice  3  compared  to 
choices  1  and  2,  but  the  effects  are  often  smaller  compared  to  the  HCRW  models.  The  HCRW 
models  further  indicate  that  increased  fragmentation  effects  should  be  expected  if  immigration 
into  MCBCL  is  reduced.  These  results  suggest  that  maintaining  migration  into  MCBCL  from 
the  Onslow  Bight  Landscape  will  help  provide  flexibility  when  determining  how  to  change  the 
allocation  of  habitat  to  relieve  training  restrictions.  Immigrants  help  reduce  influences  of  genetic 
drift  and  inbreeding  that  erode  genetic  variance  within  breeding  groups. 

Average  genetic  divergence  among  groups  (DSt)  was  more  sensitive  to  changes  in  landscape 
structure  than  Hs.  More  genetic  differences  among  groups  accrue  due  to  the  restoration  of 
dispersed  clusters  under  mitigation  choice  3.  Therefore,  the  restored  clusters  do  not  make 
equivalent  contributions  to  rates  of  migration  among  breeding  groups  as  those  that  may  be  lost 
on  the  Main  Base. 
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Figure  2.7.  LEA  credits  remaining  after  trading  on-base  clusters  for  those  on  the  Beck,  Allen, 
Lanier,  and  Greater  Sandy  Run  Area  (GSRA),  for  all  SEPMs  used  to  evaluate  changes  in  RCW 
population  dynamics.  Number  of  immigrants  per  year  from  Onslow  Bight  Landscape  were 
reduced  to  examine  how  outside  development  might  affect  habitat  trading. 
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2.4.4  Making  Cost-effective  Decisions  under  Uncertainty 

2. 4. 4.1  Decision  Making  under  Uncertainty.  Our  results  report  a  range  of 
possible  effects  that  changing  habitat  area  and  connectivity  might  have  on  RCWs  on  and 
near  MCBCL.  We  used  decision  analysis  to  more  formally  incorporate  the  role  of 
uncertainty  into  habitat  trading  or  mitigation  choices  (Raiffa  and  Schlaifer  1961;  Clemen 
and  Reilly  2001).  We  assumed  that  the  best  trade  could  be  determined  based  on  the 
number  of  LSYs  remaining  for  each  of  the  three  ecological  services.  Number  of  LSYs 
remaining  is  a  time-integrated  estimate  of  the  net  gain  or  loss  of  ecological  services  due 
to  the  change  in  landscape  structure.  The  expected  value  of  a  trade  given  uncertainty  in 
dispersal  behaviors  captured  by  the  eight  SEPMs  is  estimated  as: 

8 

EV[Trade\M]=  £  P[M]  *LSYMJrade  (eqn.  2.5) 

M-\ 

Where,  M  is  the  index  for  the  eight  uncertain  SEPMs,  P[M]  is  the  prior  probability  that 
SEPM  “M”  is  true,  LSYM,Trade  is  the  number  of  LSYs  remaining  after  the  trade  observed 
under  model  “M”.  All  eight  models  passed  through  Filter  Combination  C  from  Phase  1, 
therefore,  we  assumed  that  each  model  has  an  equal  probability  of  being  true, 

P[M]=0. 125.  The  trade  with  the  greatest  expected  value  is  referred  to  as  the  preferred 
trade  (EV[Tradep|M]).  For  PBGs  the  preferred  trade,  given  uncertainty  in  dispersal 
behaviors,  is  to  restore  all  possible  locations  outside  of  the  Main  Base  (Table  2.4).  In 
contrast,  for  both  genetic  services,  the  preferred  trade,  given  uncertainty  in  dispersal,  is  to 
restore  only  the  Lanier  tract,  which  reduced  fragmentation  effects  relative  to  the  other 
choices. 

Next  we  wish  to  know  the  expected  LSYs  remaining  assuming  that  the  best  trade  is  made 
given  the  expectations  under  each  model.  This  value  equals  the  expected  value  of 
knowing  the  true  dispersal  behaviors  (i.e.,  the  true  SEPM)  given  the  possible  trades 
available: 

8 

EV[M  | Trade]  =  £  P[M]*max(LSYM  Tracie)  (eqn.  2.6) 

M= 1 

Where,  max(LSYM, Trade)  is  the  largest  number  of  LSYs  remaining  from  the  three  possible 
trades  for  each  SEPM.  If  EV[M|Trade]  does  not  equal  EV[Trade|M]  for  all  M,  then  at 
least  one  of  the  eight  models  indicates  that  a  greater  conservation  benefit  can  be  achieved 
by  making  a  different  trade  than  suggested  when  considering  all  models  to  be  equally 
likely. 

The  expected  value  of  learning  the  true  dispersal  behaviors  can  now  be  estimated  as  the 
expected  value  of  perfect  infonnation  (EVPI  =  EV[M|Trade]  -  EV[Tradep|M]).  If  EVPI 
is  positive  then  reducing  uncertainty  may  lead  to  a  different  mitigation  choice  that 
provides  greater  conservation  benefits.  The  EVPI  is  zero  for  each  service  in  the  best  and 
worst  case  scenarios  (Table  2.4).  In  other  words,  for  genetic  services,  reducing 
uncertainty  would  likely  not  lead  to  a  decision  to  restore  other  parcels  besides  Lanier, 
because  the  expected  conservation  value  of  knowing  the  true  model  is  no  greater  than  the 
conservation  value  of  only  restoring  Lanier.  We  conclude  that  the  value  of  learning  is 
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Table  2.4.  Results  of  Decision  Analysis  applied  to  habitat  trading  using  Landscape  Equivalency  Analysis  under  uncertainty  regarding 
dispersal.  The  preferred  mitigation  choice  is  highlighted  in  gray. _ _ 


Remaining  LSYPBG 

Remaining  LSYHs 

Remaining  LSYDst 

Mitigation  Choices 

EV 

[TradeM] 

EV 

[M  Trade] 

EVPI 

EV 

[Trade  M] 

EV 

[M  Trade] 

EVPI 

EV 

[Trade  M] 

EV 

[M  Trade] 

EVPI 

Best  Case 

1.  Beck  and  Allen 

10.98 

-0.248 

-0.615 

2.  Lanier 

8.33 

-0.145 

-0.145 

0 

-0.292 

-0.292 

0 

3.  Restore  All 

13.87 

13.87 

0 

-0.385 

-0.886 

Worst  Case 

1.  Beck  and  Allen 

6.51 

-0.251 

-0.615 

2.  Lanier 

4.07 

-0.108 

-0.108 

0 

-0.241 

-0.241 

0 

3.  Restore  All 

9.46 

9.46 

0 

-0.367 

-0.875 

EV[Trade|M]  =  expected  conservation  value  of  trade  given  uncertainty  in  the  SEPMs. 
EV[M|Trade]  =  expected  conservation  value  of  knowing  the  true  SEPM  given  the  possible  trades. 
EVPI  =  EV[M| Trade]  -  EV[TradeP|M] 


Table  2.5.  Projected  management  costs  for  the  three  mitigation  choices  (Beck  and  Allen  Only,  Lanier  Only,  and  Restore  All  Possible 
outside  Main  Base).  Costs  of  land  purchase  include  contributions  from  DoD  and  Encroachment  Partners.  Management  costs  were 
taken  from  Table  16  in  ROW  Recovery  Plan  (USFWS  2003a).  _ _ _ _ 


Mitigation  Choices 

Number 

of 

Clusters 

Post  of  Land 
Acquisition 

2Initial 

Restoration 

Prescribed 

Burning 

4Cavity 

Installation 

Total  Cost 

Cost  per 
Cluster 

1.  Beck  and  Allen 

9 

$  12,097,000 

$  225,000 

$  787,500 

$  7,200 

$ 

13,116,700 

$ 

1,457,411 

2.  Lanier 

2 

$  2,032,500 

$  44,250 

$  154,875 

$  1,600 

$ 

2,233,225 

$ 

1,116,613 

GSRA5 

12 

$ 

$  300,000 

$  1,050,000 

$  9,600 

$ 

1,359,600 

$ 

113,300 

3.  Restore  All 

23 

$ 

16,709,525 

$ 

726,501 

'The  cost  of  acquiring  the  Greater  Sandy  Run  Area  was  not  included  because  this  acquisition  was  made  to  provide  land  for  military  training, 
bnitial  clearing  of  understory  at  $250/ha 

^Prescribed  burning  costs  $50/ha  to  be  applied  every  4  years  over  a  70  year  planning  horizon. 

4Assumes  that  4  cavities  are  installed  for  each  cluster  at  a  total  cost  of  $800 

5GSRA  clusters  were  not  evaluated  as  a  mitigation  choice  alone  because  they  were  the  most  isolated  clusters.  Therefore,  we  included  restoration  of  Beck,  Allen, 
and  Lanier  tracts  in  the  “Restore  All”  mitigation  choice  to  increase  connectivity  with  GSRA 
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zero  given  the  uncertainty  about  dispersal  behaviors  that  emerged  from  Phase  1  and  for 
the  trades  we  examined  in  Phase  2.  We  stress  that  this  conclusion  is  specific  to  the 
dispersal  behaviors  and  trades  that  we  considered  -  there  may  be  value  to  learning  about 
the  dispersal  uncertainties  we  considered  for  other  landscape  configurations  and  trading 
scenarios. 

2.4. 4.2  Cost-effectiveness.  LEA  organizes  habitat  trading  based  on  the  in-kind 
replacement  costs  of  the  three  ecological  services.  Replacement  costs  of  ecological 
services  were  estimated  using  the  total  cost  of  land  acquisition  and  costs  associated 
restoring  and  managing  habitat  reported  in  the  RCW  Recovery  Plan  (USFWS  2003a) 
(Table  2.5).  Therefore,  LEA  allows  us  to  incorporate  habitat  fragmentation  effects  into 
cost-effectiveness  analysis.  We  define  the  most  cost-effective  trade  as  the  one  providing 
the  greatest  conservation  benefits  at  the  lowest  cost  (Watzold  and  Schwerdtner  2005). 
Specifically,  the  most  cost  effective  trade  has  the  lowest  cost  per  LSY.  However,  the 
restoration  choices  available  were  unable  to  replace  genetic  services  and  actually  created 
a  small  net  increase  in  habitat  fragmentation.  No  trade  considered  provided  a  cost 
effective  approach  for  preventing  habitat  fragmentation.  For  discussion  purposes,  we 
present  cost  effectiveness  in  terms  of  allocating  PBGs,  an  index  of  total  habitat  area. 

Table  2.6  indicates  that  the  most  cost  effective  trade,  ignoring  habitat  fragmentation 
effects,  is  to  develop  Beck,  Allen,  Lanier,  and  GSRA.  Both  the  Letcher  model  and  the 
HCRW  model  indicate  that  this  is  the  most  cost-effective  trade  for  PBGs.  The  price  of 
credits  is  a  little  lower  under  the  Letcher  model  than  the  HCRW  models,  because  the 
Letcher  model  predicted  greater  rates  of  habitat  occupancy  in  GSRA.  The  price  of 
credits  is  higher  under  the  worse  case  development  scenarios,  because  fewer  LSYCPBG 
accrue  (Table  2.3).  Note  that  because  this  analysis  does  not  include  genetic  services,  it  is 
biased  against  mitigation  choice  2  which  resulted  in  less  fragmentation  that  the  other 
choices. 

Table  2.6.  Cost  effectiveness  of  alternative  mitigation  choices  using  the  Letcher  SEPM  and  the 

expected  value  of  the  trade  given  uncertainty  in  the  HCRW  SEPM  (EV[Trade|Mj). 


Letcher  *HCRW 

Mitigation  Choices _ LSYCFBG  Cost/LSYcPBG  LSYCPBG  Cost/LSYcPBG 


1.  Beck  and  Allen 

3.60 

Best  Case 

$3,639,726 

3.38 

$3,879,209 

2.  Lanier 

0.85 

$2,625,378 

0.73 

$3,077,291 

3.  Restore  All 

7.97 

$2,096,282 

6.28 

$2,660,246 

1.  Beck  and  Allen 

3.43 

Worst  Case 

$3,825,628 

3.16 

$4,146,555 

2.  Lanier 

0.70 

$3,195,708 

0.72 

$3,105,351 

3.  Restore  All 

7.78 

$2,147,619 

6.07 

$2,752,915 

equation  2.2 


We  found  that  the  price  of  a  LEA  credit  for  PBGs  varies  from  roughly  $2  to  $4  million 
due  to  land  values  in  coastal  North  Carolina.  Therefore,  had  decision  analysis  indicated 
that  there  was  value  in  reducing  uncertainty  regarding  dispersal  (i.e.,  positive  EVPI),  then 
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a  cost-effect  research  study  would  be  one  that  requires  less  than  (EVPI*Cost/LSY) 
dollars. 


2.5  Conclusions 

LEA  indicates  that  mitigation  choices  available  to  MCBCL  will  lead  to  different  conservation 
benefits  and  costs,  measured  as  Landscape  Service  Years  (LSYs)  for  three  ecological  services. 
The  number  of  LSYs  remaining  captured  the  contribution  individual  land  parcels  make  to 
maintaining  rates  of  habitat  occupancy,  recruitment,  and  migration  observed  under  a  clumped 
distribution  of  habitat.  If  habitat  trades  were  based  solely  on  exchanging  equal  numbers  of 
occupied  clusters,  then  MCBCL  would  be  advised  to  restore  all  possible  locations  outside  of  the 
Main  Base,  as  this  provided  the  greatest  number  of  credits  at  the  least  cost.  In  contrast,  if  habitat 
trades  were  based  on  ensuring  equivalent  rates  of  recruitment  and  migration  are  maintained  after 
a  trade,  then  none  of  the  mitigation  choices  fully  satisfy  this  condition.  However,  only  restoring 
the  Lanier  tract  led  to  the  smallest  fragmentation  effects. 

Our  analysis  suggests  that  habitat  fragmentation  effects  due  to  trading  will  likely  be  small. 

Model  assumptions  reducing  habitat  fragmentation  effects  include  implementation  of  an 
aggressive  restoration  program  on  the  Main  Base,  increasing  habitat  density,  and  constant  rates 
of  immigration  into  MCBCL  from  the  Onslow  Bight  Landscape.  These  model  conditions 
suggest  that  significant  habitat  fragmentation  effects  over  ecologically  relevant  time  scales,  such 
as  inbreeding  depression  and  increased  rates  of  local  extinction,  are  not  likely  to  result  from 
trading. 

In  general,  habitat  fragmentation  effects  were  exacerbated  by  restoring  parcels  outside  of  the 
Main  Base  (i.e.,  baseline  and  withdrawal  landscapes  were  more  functionally  equivalent  than  the 
baseline  and  mitigation  landscapes).  It  is  important  to  emphasize  that  this  was  a  consequence  of 
the  spatial  scale  of  the  analysis  and  the  fact  that  adjacency  of  clusters  increases  the  probability  of 
immigration.  The  baseline  condition  was  defined  only  based  on  the  spatial  association  of 
clusters  on  the  Main  Base.  In  the  future,  the  baseline  condition  should  be  detennined  at  a  larger 
scale  to  understand  how  habitat  trading  will  change  levels  of  fragmentation  for  the  entire  Coastal 
North  Carolina  Recovery  Unit. 

Landscape-scale  tradable  credit  systems  must  balance  uncertainty  regarding  biological  processes 
occurring  at  large  spatial  scales  with  management  costs.  Integrating  Pattern-Oriented  Modeling 
with  Decision  Analysis  provided  a  means  to  determine  if  existing  knowledge  of  behavior  and 
demography  is  sufficient  to  inform  decision  making.  Despite  existing  uncertainty  in  dispersal 
behaviors,  discussed  in  Phase  1,  it  is  unlikely  that  reducing  this  uncertainty  will  lead  to  a 
different  recommendation  for  habitat  trading,  given  the  mitigation  choices  considered. 
Uncertainty  regarding  dispersal  behaviors  did  double  the  magnitude  of  habitat  fragmentation 
effects  under  mitigation  choices  1  and  3  under  the  worst  case  development  scenario  (Figs.  2.6 
and  2.7).  Therefore,  when  restored  clusters  are  more  dispersed,  uncertainty  regarding 
attractiveness  of  land  cover  for  dispersal  may  lead  to  stronger  fragmentation  effects.  If  the 
spatial  scale  of  analysis  is  increased  to  include  the  entire  Coastal  North  Carolina  Recovery  Unit, 
then  prioritization  of  land  acquisitions  in  the  Onslow  Bight  Landscape  may  benefit  from  further 
reducing  uncertainty  regarding  dispersal. 
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3.  Phase  3:  Exploring  Habitat-Fitness  Relationships  with  Pre-Existing  Data 

3.1  Objective 

The  objective  of  this  phase  is  to  derive  a  new  recruitment  submodel  that  can  be  initialized 
without  knowing  the  age-structure  of  the  population,  so  that  the  SEPM  can  be  applied  at  a 
broader  spatial  scale.  We  describe  how  habitat-fitness  relationships  can  be  determined  using  the 
foraging  habitat  data  collected  in  1998-2000  by  J.H.  Carter  &  Assoc,  and  MCBCL  (GEC  2000). 
These  data  do  not  meet  all  the  requirements  of  USFWS’s  updated  Foraging  Habitat  Guidelines 
(USFWS  2003a),  because  they  did  not  collect  infonnation  on  groundcover.  While  it  is  important 
to  use  USFWS  guidelines  for  managing  RCWs,  new  methods  that  lead  to  cost-effective 
approaches  of  landscape-scale  management  should  continue  to  be  investigated.  We  introduce  a 
new  approach  that  summarizes  forest  attributes  within  different  lag  distances  from  the  cavity 
trees.  Our  analysis  also  includes  the  effect  of  year  into  habitat- fitness  relationships.  Year  of 
observation  was  included  because  there  seems  to  be  large  yearly  fluxes  in  recruitment  that  may 
be  due  to  climatic  effects  such  as  rainfall. 

3.2  Background 

A  major  benefit  of  using  individual-based,  spatially-explicit  population  models  (SEPMs)  is  their 
ability  to  simulate  interactions  among  population  processes  and  changes  in  landscape  pattern  that 
result  from  habitat  trading.  We  might  call  such  models  mechanistic  and  non-equilibrium  models. 
If  we  wish  to  trade  habitat  across  the  entire  geographic  extent  of  the  Coastal  North  Carolina 
population,  then  the  model’s  ability  to  predict  population  dynamics  across  this  entire  area  should 
be  incorporated  into  decision  making.  In  other  words,  population  monitoring  data  collected 
outside  of  Marine  Corps  Base  Camp  Lejeune  (MCBCL)  should  serve  as  patterns  used  to 
initialize  the  model  as  well  as  to  estimate  and  reduce  uncertainty,  as  described  in  Phases  1  and  2. 
Model  calibration  and  validation  should  be  an  integral  part  of  any  decision  making  process  that 
relies  on  SEPMs,  due  to  their  reliance  on  uncertain  parameters  and  structural  assumptions 
(Beisinger  and  Westphal  1998;  Wiegand  et  al.  2004). 

As  previously  presented,  the  SEPMs  used  in  this  study  and  by  DoD  and  DOI  for  managing  Red- 
Cockaded  Woodpeckers  (RCWs)  rely  on  information  regarding  the  age  structure  and  social 
structure  of  RCW  populations  to  initialize  the  model  and  simulate  population  dynamics.  While 
this  quality  of  data  is  available  for  MCBCL,  it  is  not  available  for  other  management  units  of  the 
Coastal  North  Carolina  population.  For  example,  data  on  the  age  structure  of  the  entire 
population  is  not  available.  The  RCW  SEPMs  (Phases  1  and  2)  predict  the  number  of  fledglings 
in  a  territory  based  on  the  ages  of  the  male  and  female  breeder  and  number  of  helpers  (Letcher  et 
al.  1998).  The  age-  and  helper-based  recruitment  model  is  based  on  a  statistically  significant 
regression  equation  (p<0.0001),  but  only  accounts  for  8%  of  the  observed  variance  in  fledglings 
(D.  Bruggeman,  Phase  1). 

It  may  be  more  cost-effective  to  collect  data  on  habitat  attributes  than  to  continue  bird  banding 
studies  that  lead  to  estimates  of  age  structure,  especially  as  the  population  size  increases  toward 
recovery.  Previous  attempts  to  establish  a  statistically  significant  relationship  between  the 
number  of  fledglings  produced  and  vegetative  characteristics  of  Red-Cockaded  Woodpecker 
(RCW)  territories  (i.e.,  fitness-habitat  relations)  have  not  been  successful  (James  et  al.  2001; 
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Walters  et  al.  2002).  However,  many  biologists  have  observed  that  groups  in  mature  Longleaf 
pine  forests  that  are  regularly  subjected  to  burning  tend  to  have  more  productive  groups. 

Convery  (2002)  did  successfully  relate  recruitment  to  habitat  attributes  in  23  territories  on 
MCBCL.  Based  on  the  observations  of  James  et  al.  (1997;  2001),  Convery  (2002)  supplemented 
previous  surveys  of  RCW  foraging  habitat  (GEC  2000)  by  collecting  data  on  ground  cover, 
because  it  is  strong  indicator  of  fire  history.  Convery  (2002)  also  conducted  behavioral 
observations  to  estimate  size  and  location  of  home  ranges.  With  this  approach,  Convery  derived 
statistically  significant  relationships  between  forest  attributes,  based  on  USFWS’s  updated 
Foraging  Habitat  Guidelines  (USFWS  2003a),  and  recruitment  using  linear  regression.  The 
work  of  Convery  provides  hope,  but  the  efforts  required  to  estimate  home  ranges  and  quantify 
ground  cover  have  not  been  replicated  across  MCBCF  and  would  be  time  consuming  to 
implement  across  the  entire  Coastal  North  Carolina  population. 

3.3  Methods 


3.3.1  Foraging  Habitat  Data  for  MCBCL 

Concerns  regarding  changes  in  forest  structure  on  MCBCL  due  to  hurricanes  and  pine  beetle 
infestations  led  to  a  major  effort  to  characterize  RCW  foraging  habitat  across  the  entire  base 
from  1998-2000  (GEC  2000).  Initially,  5,670  0.04-hectare  plots  were  laid  out  across  MCBCL  on 
a  standard  grid  of  five  chains  north  to  south  by  12  chains  east  to  west.  Sampling  focused  on 
managed  pine  and  pine -hardwood  stands  that  were  at  least  30  years  old,  but  some  younger  stands 
were  also  sampled.  A  total  of  5,374  forested  plots  were  characterized. 

During  the  foraging  habitat  survey  the  number  of  stems  per  hectare  was  estimated  in  5  cm  dbh 
classes,  for  pines  and  hardwoods.  Following  Convery  (2002),  we  condensed  these  data  into 
three  size  classes  for  pines  to  parallel  USFWS  Foraging  Habitat  Guidelines  (Table  3.1). 

Similarly,  we  condensed  data  for  hardwoods  into  two  size  classes. 

3.3.2  Spatially  Integrating  Foraging  Habitat  Data 

Home  ranges  for  all  territories  on  MCBCL  have  not  been  estimated.  To  characterize  habitat 
attributes  corresponding  to  individual  RCW  groups,  we  developed  a  new  approach.  USFWS 
Foraging  Habitat  Guidelines  specifies  that  foraging  habitat  should  be  within  800  m  of  a  cluster’s 
center.  Convery  (2002)  estimated  average  home  range  size  as  80.2  ha  based  on  23  territories  on 
MCBCL  by  analyzing  the  movement  data  using  the  Animal  Movement  Analysis  extension  in 
ArcView  GIS.  This  suggests  that  foraging  habitat  tends  to  be  within  800  m  even  if  the  cavity 
trees  were  located  at  a  border  of  a  home  range.  We  summarized  foraging  habitat  data  using  a 
moving  window  approach  to  integrate  forest  attributes  at  different  lag  distances  from  cluster 
centers.  This  approach  helps  control  uncertainty  regarding  the  borders  of  home  ranges  and  the 
influence  of  forest  fragmentation  on  habitat  quality.  It  may  lead  to  a  useful  approach  for 
determining  how  loss  of  forest  near  clusters  may  affect  population  dynamics.  Three  sizes  of 
windows  were  employed:  100  ha,  400  ha,  and  1,225  ha  (Table  3.1).  These  windows  correspond 
to  500  m,  1,000  m,  and  3,500  m  lag  distances  from  the  cluster  center.  Note,  the  last  lag  distance 
corresponds  to  estimated  3.5  km  search  radius  of  helpers  (Pasinelli  and  Walters  2002).  To 
characterize  habitat  attributes  within  a  window,  we  calculated  the  average  density  of  each  size 
class  of  pines  and  hardwoods: 
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1  5 

PSw  =  —  X  PineSmall(stems / ha)  j  (eqn.  3.1) 

W7=l 

where  PS  means  small  pines  (Table  3. 1),  w  is  the  size  of  the  window  in  hectares,  and  s  is  the 
number  of  times  small  pines  were  detected  within  a  window. 

3.3.3  Regression  Approach 

Simple  linear  regression  was  used  to  parallel  other  studies  (Convery  2002;  James  et  al.  2001; 
Walters  et  al.  2002).  Figure  3.1  indicates  that  the  number  of  fledglings  per  cluster  approximate  a 
nonnal  distribution,  so  simple  linear  regression  is  appropriate.  However,  unlike  previous  studies 
we  chose  to  include  the  intercept  tenn  during  model  fitting.  Biologically  it  makes  sense  to 
exclude  the  intercept  because  areas  with  no  trees  will  not  contain  fledglings.  However,  omitting 
the  intercept  inflates  the  estimate  of  model  fit,  R  ,  and  increases  prediction  bias.  Ideally,  the 
relationship  between  group  fitness  and  the  independent  variables  predicted  by  the  model  should 
remove  the  intercept  from  the  model  (i.e.,  the  intercept  should  not  be  significantly  different  from 
zero). 
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Figure  3.1.  Frequency  distribution  of  fledglings  per  group  (fitness)  observed  from  1998  to  2000 
on  MCBCL. 

Independent  variables  used  in  the  regression  were  based  on  the  foraging  habitat  survey  described 
previously  and  other  variables  that  may  be  easily  collected  across  the  extent  of  the  Coastal  North 
Carolina  population.  It  is  well  accepted  that  clusters  with  more  helpers  tend  to  produce  more 
fledglings  (Walters  et  al.  2002;  James  et  al.  2001;  Letcher  et  al.  1998),  so  we  supplemented 
habitat  variables  with  observed  number  of  male  and  female  helpers  per  cluster.  Also,  Charles 
Clarkson  (VPI&SU),  who  manages  the  RCW  monitoring  program  at  MCBCL,  indicated  that 
rainfall  patterns  likely  have  a  large  impact  on  recruitment.  To  conduct  a  preliminary  analysis  of 
this  hypothesis,  we  included  year  of  observation  as  an  independent  variable  (Table  3.1). 
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Table  3.1.  Derivation  of  independent  variables  for  use  in  habitat-fitness  relationships. 


'Abbreviation 

Independent 

Variable 

Converting  data  from 
MCBCL  Foraging 
Habitat  Survey 

Corresponding  USFWS 
(2003)  Recommendation 

PL1,  PL4, 

PL12 

Large  Pines 

No.  pine  stems/ha  >  35 
cm  dbh 

2a.  >  45  stems/ha  of  pine  at 
least  35  cm  dbh 

PM1,  PM4, 

PM  12 

Medium  Pines 

No.  pine  stems/ha  25.4  to 
35  cm  dbh 

2b.  basal  area  of  pines  25.4  to 

35  cm  dbh  between  0  and  9.2 
m  /ha 

PS1.PS4, 

PS12 

Small  Pines 

No.  pine  stems/ha  <  25.4 
cm  dbh 

NA 

HL1,  HL4, 

Large 

No.  hardwood  stems/ha  > 

2g.  Canopy  hardwoods  <10% 

HL12 

Hardwoods 

25.4  cm  dbh 

of  canopy  in  longleaf  stands, 
and  <  30%  of  canopy  in 
loblolly  stands 

HS1,  HS4, 

Small 

No.  hardwood  stems/ha  < 

2f.  Hardwood  midstory  <  2.1 

HS12 

Hardwoods 

25.4  cm  dbh 

m 

HP 

No.  Helpers 

NA 

NA 

YR 

Year  of 
Observation 

NA 

NA 

e.g.,  PL1  =  density  of  large  pines  within  100-ha  window. 


Statistical  analyses  were  completed  in  the  software  R  version  2.2.0  (Venables  and  Smith  2005). 
Three  sets  of  full  regression  equations  were  compiled.  The  first  included  only  the  habitat 
variables.  The  second  included  all  the  habitat  variables  and  number  of  helpers  per  group.  The 
third  included  habitat,  helpers,  and  year  of  observation. 

Stepwise  model  selection  was  performed  using  the  Akaike  Information  Criterion  (AIC).  AIC 
uses  the  Kullback-Leibler  distance  to  quantify  the  discrepancy  between  the  distribution 
describing  the  true  outcomes  of  an  environmental  process  and  the  distribution  of  outcomes 
predicted  by  the  model  (Burnham  and  Anderson  2002).  Because  the  distribution  of  true 
outcomes  can  never  be  known,  the  AIC  method  estimates  the  Kullback-Leibler  distance  by  using 
maximum  likelihood  to  fit  observed  data  to  the  model.  Unlike  the  more  traditional  statistical 
estimate  “R  ”,  which  only  estimates  model  fit  based  on  the  proportion  of  the  sum  of  squares 
explained  by  a  model,  AIC  incorporates  the  number  of  parameters  into  the  assessment  of  model 
fit  (Richards  2005).  The  most  parsimonious  model,  capable  of  explaining  the  largest  amount  of 
variance  in  the  data  with  the  fewest  parameters,  is  given  the  lowest  AIC  score.  The  “stepAIC” 
procedure  in  the  R  2.2.0  software  package  was  used  to  determine  the  model  with  the  lowest  AIC 
score  using  forward  and  backward  elimination.  However,  when  models  differ  in  AIC  scores  by 
2  units  or  less,  then  we  say  that  the  two  models  are  equivalent  (Richards  2005).  Therefore,  we 
also  report  the  most  parsimonious  model  within  2  AIC  units  of  the  model  selected  by  the 
stepAIC  procedure. 
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3.4  Results 

Regressing  group  fitness  against  only  habitat  variables  did  not  lead  to  a  statistically  significant 
relationship  (p=0.063;  Model  A,  Table  3.2).  However,  AIC  model  selection  indicated  that  four 
habitat  variables  explained  3.3%  of  the  variation  in  group  fitness.  Density  of  small  pines  within 
a  100-ha  window  about  the  cluster  center  was  positively  related  to  group  fitness.  In  contrast,  the 
density  of  medium  pines  within  this  same  window  was  negatively  related  to  group  fitness. 

Within  a  400-ha  window,  the  density  of  small  pines  was  negatively  related  to  group  fitness. 
Similarly,  within  a  1,225-ha  window  the  density  of  small  hardwoods  is  negatively  related  to 
group  fitness.  Even  though  the  model  is  not  significant,  these  results  suggest  that  group 
productivity  may  be  sensitive  to  forest  structure  well  outside  of  home  ranges.  When  we 
simplified  the  model  further  without  increasing  the  AIC  score  greater  than  2  units,  we  found  that 
only  density  of  medium  pines  within  100  ha  was  retained. 

Model  B  summarizes  the  results  of  including  the  number  of  helpers  and  stem  densities  for  pines 
and  hardwoods  at  different  lag  distances  (Table  3.2).  Model  B  explained  8.9%  of  the  variance  in 
group  fitness  and  resulted  in  a  large  decrease  in  the  AIC  score  compared  to  Model  A.  AIC 
model  selection  resulted  in  excluding  all  habitat  variables  except  the  density  of  large  hardwoods 
within  a  100-ha  window  and  the  density  of  small  pines  and  medium  pines  within  a  1,225 -ha 
window  (Model  B,  Table  3.2).  The  model  predicts  that  density  of  large  hardwoods  near  the 
cavity  trees  is  negatively  related  to  group  fitness.  At  the  largest  window  size,  the  model 
indicates  that  the  density  of  small  pines  is  positively  related  with  group  fitness,  but  density  of 
medium  pines  is  negatively  related  to  group  size.  The  correlation  between  PS  12  and  PM12  was 
-0.62,  which  is  sufficiently  small  to  include  both  variables.  Further  simplifying  Model  B 
indicated  that  an  equivalent  model  would  be  one  that  includes  only  number  of  helpers. 

When  the  effect  of  year  was  added  into  AIC  model  selection,  we  observed  a  significant 
improvement  in  model  fit  and  the  same  habitat  variables  were  retained  (Model  C,  Table  3.2). 
Helpers  still  play  a  dominant  role  in  determining  the  number  of  fledglings,  but  the  influence  of 
year  was  more  significant  than  the  habitat  variables.  With  the  inclusion  of  more  independent 
variables  the  adjusted  R2  value  increased  from  0.089  (Model  B)  to  0. 12  in  Model  C,  which  is 
expected  because  we  have  more  variables  present  to  account  for  variance.  However,  the  AIC 
score  dropped  by  5  units,  indicating  that  the  increased  model  complexity  does  provide  useful 
information  regarding  the  relationship  between  environmental  variables  and  fitness.  However, 
the  predicted  distribution  of  fledglings  per  group  (Fig.  3.2)  underestimates  the  observed  number 
of  fledglings  per  group  (Fig.  3.1).  When  we  further  simplified  the  model,  only  number  of 
helpers  and  year  was  retained. 
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Table  3.2.  Habitat- fitness  relationships  predicted  by  linear  regression  and  AIC  model  selection. 
Data  analyzed  represent  154  Potential  Breeding  Groups  observed  from  1998  to  2000. 


^odel  A 

p-value 

Ad  j  R2  AIC 

Fit  =  2.50  +  0.0079  PS1  -  0.036  PM1  -  0.0102  PS4  - 

-  0.019  HS12 

0.063 

0.033  426 

Parameter 

p-value 

t-value 

Intercept 

<0.0001 

8.23 

PS1 

0.125 

1.54 

PM1 

0.038 

-2.10 

PS4 

0.107 

-1.62 

HS12 

0.063 

-1.88 

Model  A  -  Most  Parsimonious 

p-value 

Adj  R2  AIC 

Fit  =  2.02  -  0.03  PM1 

0.053 

0.018  427 

^odel  B 

p-value 

Adj  R2  AIC 

Fit  =  1.93  +  0.247  HP  -  0.238  HL1  +  0.0085  PS  12  - 

0.043  PM  12 

0.0014 

0.089  418 

Parameter 

p-value 

t-value 

Intercept 

<0.0001 

4.66 

HP 

0.003 

2.99 

HL1 

0.042 

-2.06 

PS12 

0.163 

1.40 

PM12 

0.043 

-2.04 

Model  B  -  Most  Parsimonious 

p-value 

Adj  R2  AIC 

Fit  =  1.40  +  0.28  HP 

0.0007 

0.067  419 

^odel  C 

p-value 

Adj  R2  AIC 

Fit=1.51+0.252HP+0.228YR-0.253HL1+0.009PS12 

-0.047PM12 

0.0002 

0.120  413 

Parameter 

p-value 

t-value 

Intercept 

0.0007 

3.44 

HP 

0.0023 

3.11 

YR 

0.013 

2.52 

HL1 

0.028 

-2.22 

PS12 

0.142 

1.48 

PM12 

0.028 

-2.22 

Model  C  -  Most  Parsimonious 

p-value 

Adj  R2  AIC 

Fit  =  0.96  +  0.29  HP  +  0.21  YR 

0.0002 

0.093  415 

’inclusion  of  all  other  independent  variables  only  increased  AIC  scores. 
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Predicted(Fitness)  Residual(Fitness) 

Figure  3.2.  (A)  Frequency  distribution  of  fledglings  per  group  (fitness)  predicted  by  Model  C. 
(B)  Frequency  distribution  of  residuals  from  Model  C. 


3.5  Discussion 

We  have  demonstrated  that  integrating  data  describing  vegetative  characteristics  at  different 
spatial  scales,  number  of  helpers,  and  differences  among  years  can  lead  to  a  better  prediction  of 
fledglings  than  the  age-  and  helper-based  recruitment  model  already  used  for  decision  making. 
However,  there  is  still  much  room  for  improvement,  because  at  most  only  12%  of  the  variation  in 
recruitment  was  explained  by  the  model.  Two  observations  made  by  the  current  work  suggest 
that  further  research  might  improve  our  prediction  of  fledglings.  First,  density  of  vegetation 
outside  of  home  ranges  is  important.  Second,  yearly  variation  in  environmental  conditions 
should  be  taken  into  account.  These  are  further  discussed  below. 


The  importance  of  continuous  old  growth  Longleaf  pine  forests  has  been  well  articulated  by 
other  researchers  (Conner  and  Rudolph  1991;  Cox  et  al.  2001;  James  et  al.  2001;  Walters  et  al. 
2002).  However,  quantitative  analysis  of  forest  fragmentation  effects  have  been  limited  to  an 
800-m  lag  distance  from  the  cluster  center  (Conner  and  Rudolph  1991;  Cox  et  al.  2001;  James  et 
al.  2001),  based  on  USFS  and  USFWS  guidelines.  In  contrast,  this  study  included  a  500-m, 
1,000-m,  and  3, 500-m  lag  distance,  corresponding  to  100-ha,  400-ha,  and  1,225-ha  square 
window.  This  approach  provided  interesting  results.  Within  the  100-ha  and  1,225-ha  windows, 
AIC  model  selection  suggested  that  the  density  of  small  pine  was  positively  related  to  fitness. 
However,  high  densities  of  small  pines  in  a  400-ha  window  was  negatively  related  to  fitness. 
While  the  coefficients  describing  the  influence  of  small  pines  were  not  statistically  significant 
(p>0.05),  AIC  model  selection  indicated  that  they  do  contribute  to  model  fit. 

James  et  al.  (2001)  found  that  density  of  small  pines  (15-25  cm  dbh)  was  negatively  related  to 
number  of  fledglings,  but  that  the  difference  in  large  to  small  pines  was  positively  related  to 
fitness.  Walters  et  al.  (2002)  suggests  that  RCWs  tend  to  avoid  areas  with  small  pines,  are 
indifferent  to  medium  pines,  and  prefer  large  pines  during  foraging.  These  results  seem  to 
contradict  our  observations.  We  attempt  to  reconcile  these  observations  below,  by  stressing  the 
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conclusion  of  James  et  al.  (2001)  that  the  structure  of  pine  forest  is  critical,  not  just  the  densities 
of  individual  sizes  classes  of  pines. 

Recall  that  almost  all  plots  characterized  on  MCBCL  were  managed  pine  stands  30  years  or 
older.  Therefore,  presence  of  small  pines  near  cavity  trees  (i.e.,  within  a  100-ha  window)  is 
likely  indicative  of  natural  pine  regeneration  and  fire  history.  Management  of  longleaf  pine 
forests  through  prescribed  burning  leads  to  uneven-aged  forest  stands,  in  which  low  densities  of 
larger  pines  are  interwoven  with  small  areas  with  higher  densities  of  smaller  pines  able  to 
regenerate  in  gaps  created  by  the  loss  of  an  older  tree,  absent  competition  from  hardwoods  in  the 
understory  (James  et  al.  2001).  In  contrast,  high  density  of  small  pines  within  a  400-ha  window, 
farther  from  cavity  trees,  would  suggest  the  presence  of  even  aged  forest  due  to  a  recent  planting. 
Similarly,  high  density  of  medium  pines  within  a  100-ha  window  may  be  indicative  of  an  even 
aged  forest  near  cavity  trees,  which  the  model  predicts  will  lead  to  lower  fitness.  Model  A 
suggests  that  density  of  small  hardwoods  within  1,225-ha  was  negatively  related  to  fitness.  This 
suggests  that  maintaining  historical  fire  regimes  well  outside  of  a  group’s  home  range  could 
contribute  to  group  fitness. 

When  number  of  helpers  were  included,  variables  describing  pine  in  the  smaller  two  windows 
dropped  out  of  the  model.  This  suggests  that  number  of  helpers  is  a  good  surrogate  for  attributes 
of  pine  structure  describing  good  quality  habitat  close  to  cavity  trees  and  home  ranges.  Model  B 
does  include  density  of  larger  hardwoods  within  a  100-ha  window.  The  negative  relationship 
between  fitness  and  HL1  corresponds  to  observations  of  other  researchers  that  RCWs  avoid 
forested  areas  with  large  hardwoods  when  foraging  (Walters  et  al.  2002).  The  positive 
relationship  between  density  of  small  pines  in  the  largest  window  (1,225-ha)  and  fitness  in 
Model  B  may  be  indicative  natural  pine  regeneration  at  large  spatial  scales,  but  this  variable  does 
not  contribute  greatly  to  model  fit  (t-vale  =  1.4).  The  negative  relationship  between  density  of 
medium  pine  in  the  largest  window  and  fitness  may  reflect  the  presence  of  large  even  aged 
stands.  It  is  curious  that  these  larger  scale  habitat  variables  are  introduced  only  after  including 
the  effect  of  helpers.  Perhaps  the  structure  of  pine  forests  within  a  1,225-ha  window  may  affect 
the  probability  that  a  breeding  vacancy  is  filled  by  a  floater?  Recall  that  Phase  1  attempted  to 
characterize  habitat  effects  on  dispersal,  but  used  much  cruder  estimates  of  forest  structure  (i.e., 
size  class  distributions  were  ignored). 

Yearly  environmental  variation  contributed  considerably  to  fitness  (Model  C).  The  t-value  for 
year  provided  the  third  largest  absolute  value.  Further,  inclusion  of  the  year  effect  helped  to 
increase  the  significance  of  all  other  independent  variables,  suggesting  a  possible  interaction 
between  year  and  habitat  quality.  To  the  best  of  our  knowledge  yearly  variation  in  recruitment 
has  not  been  considered  in  previous  studies.  When  multi-year  data  were  available,  researchers 
averaged  fitness  within  territories  across  years  before  perfonning  regression  analysis  (Convery 
2002;  Walters  et  al.  2002).  Arguably,  including  year  as  a  factor  is  not  the  best  surrogate  for 
environmental  variables  that  could  affect  recruitment  such  as  rainfall  and  temperature.  One 
confounding  factor  could  be  that  the  number  of  territories  increased  from  46  in  1998  to  56  in 
2000.  However,  the  presence  of  these  new  clusters  could  have  reduced  the  number  of  helpers 
per  group  due  to  new  breeding  vacancies.  We  strongly  encourage  future  research  that  examines 
how  climatic  variables  might  contribute  to  recruitment. 
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In  the  SERDP  SEED  grant  contract  modified  proposal,  we  indicated  that  the  new  habitat-fitness 
relationship  derived  here  would  be  included  into  the  RCW  SEPM.  However,  this  investigatory 
research  discovered  that  yearly  variation  in  environmental  factors  is  key  for  providing  a 
relationship  that  is  at  least  as  good  as  the  age-  and  helper-based  recruitment  model  currently 
being  used.  Due  to  uncertainty  regarding  what  those  factors  are  (e.g.,  rainfall  levels  during 
which  months),  we  were  unable  to  update  the  SEPM. 

3.6  Conclusion 

Linear  regression  modeling  including  a  moving  window  approach  and  yearly  variation  in 
recruitment  did  improve  the  prediction  of  group  fitness.  Observations  including  numbers  of 
helpers  within  groups,  forest  structure  at  multiple  scales,  and  climatic  variables  could  be 
collected  without  instituting  a  bird  banding  study  across  the  Coastal  North  Carolina  population 
ofRCWs. 

The  forest  attributes  included  in  this  study  were  collected  at  great  expense  to  MCBCL  and  were 
not  being  used  for  management.  Replicating  this  effort  in  Croatan  National  Forest,  Holly  Shelter 
State  Game  Lands,  and  smaller  areas  considered  for  acquisition  under  the  Encroachment 
Partnering  program  may  be  expensive.  Remotely  sensed  data  describing  vegetative  structure 
may  also  be  integrated  with  climatic  data  and  may  be  more  cost-effective.  In  Phase  3,  we  have 
been  able  to  devise  an  approach  to  bring  value  to  existing  data  while  directing  future  data 
collection  efforts. 
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4.  Overall  Conclusions  and  Recommendations 

4.1  Summary  of  Significant  Results 

The  use  of  Encroachment  Partnering  to  relieve  training  restrictions  due  to  the  presence  of 
endangered  species  highlights  the  importance  of  understanding  how  changing  landscape  patterns 
may  affect  recruitment  and  migration.  Individually-based,  spatially-explicit  population  models 
(SEPMs)  can  be  used  to  determine  how  changing  habitat  area  and  /  or  connectivity  will  affect 
population  viability.  These  models  allow  interaction  between  rates  of  recruitment  within  patches 
with  rates  of  migration  among  patches.  We  applied  new  techniques  to  improve  how  SEPMs  can 
be  used  for  Encroachment  Partnering.  The  major  methodological  advances  and  conclusions  are 
summarized  below: 

1 .  Pattern-Oriented  Modeling  (POM)  provides  a  strong  inference  approach  for  reducing 
uncertainty  regarding  parameter  values  in  SEPMs.  Out  of  100,000  alternative  dispersal 
submodels  tested,  only  8  were  able  to  simultaneously  reproduce  all  ROW  patterns 
observed  on  MCBCL. 

2.  Of  the  8  parameter  sets  remaining,  estimates  of  landscape  friction  and  turning  angles 
varied  greatly.  We  proposed  two  possible  reasons  for  this  result.  One,  POM  was  applied 
over  too  small  an  area.  In  other  words,  the  ability  of  simulation  model  to  reproduce 
patterns  observed  on  MCBCL  was  not  affected  by  large  differences  in  land  use.  Two,  the 
model  does  not  accurately  describe  how  dispersing  birds  make  decisions  to  maximize 
survivability  and  to  find  a  breeding  vacancy  (i.e.,  structural  uncertainty  remains  in  the 
model). 

3.  Landscape  Equivalency  Analysis  (LEA)  provides  a  quantitative  approach  for  determining 
how  much  regulatory  relief  can  be  achieved  from  specific  Encroachment  Partnering 
decisions.  By  calibrating  the  new  dispersal  submodel  to  patterns  observed  on  MCBCL, 
we  found  that  Encroachment  Partnering  may  lead  to  less  regulatory  relief  than  would  be 
projected  using  the  dispersal  submodel  found  in  the  Letcher  model. 

4.  Integrating  uncertainty  analysis  (POM)  with  a  landscape-scale  tradable  credit  system 
(LEA)  indicated  that  we  have  sufficient  knowledge  to  prioritize  habitat  restoration 
decisions. 

5.  We  used  decision  analysis  to  describe  how  DoD  can  determine  when  it  is  cost  effective  to 
forego  cooperative  purchase  of  off-base  parcels  and  collect  more  biological  data  to  ensure 
that  Encroachment  Partnering  will  lead  to  regulatory  relief  without  decreasing  population 
viability. 

6.  Habitat-fitness  relationships  were  found  by  applying  a  moving  window  approach  to 
quantify  forest  structure  at  multiple  spatial  scales.  We  found  that  yearly  variation  in 
environmental  conditions  contribute  significantly  to  observed  variation  in  recruitment. 
This  approach  accounted  for  more  of  the  total  variance  in  recruitment  (Adj.  R  =0.12) 
than  the  age-base  recruitment  submodel  (Adj.  R  =0.08)  currently  used  by  the  SEPM. 

In  conclusion,  by  examining  habitat  trading  opportunities  on  and  near  MCBCL,  we  found  that 
POM  reduced  enough  uncertainty  regarding  dispersal  to  prioritize  trading  decisions.  However, 
we  also  observed  that  examining  habitat  trades  at  this  spatial  extent  always  resulted  in  increasing 
habitat  fragmentation  effects. 
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It  is  critical  to  determine  if  Encroachment  Partnering  will  increase  or  decrease  habitat 
fragmentation  for  the  entire  recovery  population.  MCBCL  is  one  of  three  Management  Units 
(MUs)  within  the  Coastal  North  Carolina  RCW  recovery  population.  The  other  two  include 
Croatan  National  Forest  and  Holly  Shelter  Game  Lands,  and  the  three  MUs  are  found  in  the 
Onslow  Bight  Landscape.  Habitat  trades  that  increase  habitat  fragmentation  relative  to  the 
current  population  on  MCBCL  may  actually  decrease  fragmentation  across  the  entire  Coastal 
North  Carolina  recovery  population  inhabiting  the  Onslow  Bight  landscape.  Further,  the 
remaining  uncertainty  in  dispersal  would  likely  affect  habitat  trading  decisions  over  the  entire 
Onslow  Bight  Landscape.  Given  the  high  rate  of  sprawl  in  coastal  North  Carolina,  MCBCL  has 
initiated  an  aggressive  Encroachment  Partnering  program.  For  instance,  MCBCL  is  working  to 
protect  lands  for  RCWs  that  are  closer  to  Holly  Shelter  Game  Lands  than  MCBCL’s  Main  Base. 
Determining  how  much  regulatory  relief  can  be  achieved  through  such  actions  indicates  that 
increasing  the  spatial  scale  of  our  model  is  a  critical  area  of  future  research. 

4.2  Priorities  for  Future  Research 

4.2.1  Habitat  Trading  for  RCWs  in  Onslow  Bight  (Bruggeman,  JR  Walters  et  al.) 

In  order  to  detennine  if  DoD  can  spread  RCW  liability  across  the  Onslow  Bight  landscape 
without  jeopardizing  recovery,  we  need  to  understand  how  landscape  patterns  affect  recruitment 
and  migration  at  a  broad  geographic  extent.  Results  from  Phase  3  (habitat-fitness  analysis)  are 
currently  being  explored  in  collaboration  between  Douglas  Bruggeman  and  Jeffrey  Walters 
under  SI- 1472.  Climatological  data  have  been  collected  for  Fort  Bragg  and  MCBCL  and  will  be 
integrated  with  the  foraging  habitat  quality  data  Dr.  Walters’  group  has  been  compiling.  This 
effort  should  shed  light  on  the  cause  of  the  strong  “year-effect”  observed  in  Phase  3.  Dispersal 
uncertainty  is  being  addressed  in  a  radiotelemetry  study  (SI- 1471)  through  collaboration  between 
Aaron  Moody  and  Jeffrey  Walters.  Douglas  Bruggeman  has  recently  shared  results  of  the  POM 
analysis  of  dispersal  uncertainty  (Phase  2)  to  help  guide  the  integration  of  the  radiotelemetry 
studies  with  the  RCW  SEPM.  Despite  this  increased  collaboration,  there  is  still  one  critical 
research  question  not  being  addressed  under  any  existing  proposal: 

What  are  the  historic  and  current  rates  of  migration  among  the  three  Management  Units  in 
the  Onslow  Bight  (Croatan  National  Forest,  MCBCL,  and  Holly  Shelter  Game  Lands)? 
Understanding  historic  and  current  landscape  function  at  this  scale  is  critical  for 
determining  how  much  regulatory  relief  can  be  achieved  through  Encroachment 
Partnering. 

In  order  capture  habitat  fragmentation  effects,  LEA  assigns  a  conservation  value  to  trades  at  a 
local  scale  based  on  marginal  changes  in  ecological  services  at  the  landscape  (regional)  scale, 
relative  to  service  levels  expected  in  an  unfragmented  or  baseline  landscape.  We  defined 
baseline  for  MCBCL  as  the  apportionment  of  genetic  variance  within  and  among  breeding 
groups  afforded  by  their  recovery  plan.  By  applying  the  concept  to  the  entire  Onslow  Bight,  we 
might  define  baseline  as  the  apportionment  of  genetic  variance  within  and  among  breeding 
groups,  averaged  across  all  three  management  units  (MUs),  afforded  by  reaching  property- 
specific  recovery  goals  (USFWS  2003a).  In  contrast  to  the  baseline  defined  for  MCBCL,  the 
baseline  for  the  Onslow  Bight  does  not  provide  the  most  clumped  distribution  of  habitat 
possible,  merely  the  most  clumped  distribution  possible  on  public  lands.  We  propose 
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research  to  determine  if  this  provides  adequate  levels  of  connectivity  for  RCWs.  The  work 
would  also  determine  if  restoration  of  private  land  would  increase  population  viability,  justifying 
the  use  of  Encroachment  Partnering  to  provide  regulatory  relief. 

Landscape  Genetic  Study  of  Coastal  North  Carolina  RCW  Recovery  Population 

We  propose  to  conduct  a  landscape  genetics  study  in  the  Onslow  Bight  landscape.  Tissues  will 
be  collected  from  RCWs  located  on  MCBCL,  Croatan  National  Forest,  and  Holly  Shelter.  These 
data  will  be  used  to  address  the  following  objectives. 

1 .  Estimate  levels  of  connectivity  prior  to  habitat  loss  and  fragmentation.  The  tissues 
will  be  assayed  for  genetic  variability  using  allozymes  and  micro  satellites,  because  these 
two  classes  of  nuclear  DNA  have  very  different  properties.  Microsatellites  mutate 
frequently  but  are  selectively  neutral.  In  contrast,  allozymes  mutate  much  less  frequently 
and  the  genotypes  may  be  associated  with  differences  in  fitness  among  individuals, 
breeding  groups,  or  MUs.  Therefore,  we  can  estimate  habitat  fragmentation  effects  at 
different  time  scales  by  contrasting  the  spatial  apportionment  of  allozymes  to 
microsatellites. 

The  spatial  apportionment  of  microsatellite  markers  provides  an  opportunity  to  test  the 
implicit  assumption  of  the  RCW  recovery  plan  that  land  use  changes  in  between  MUs 
will  not  affect  recovery.  In  contrast,  the  spatial  apportionment  of  allozyme  markers  will 
be  descriptive  of  migration  rates  prior  to  fragmentation.  We  predict  that  the  baseline 
levels  of  genetic  services  provided  by  the  RCW  recovery  plan  will  likely  be  an 
inadequate  representation  of  an  unfragmented  condition.  Knowledge  of  historic 
apportionment  of  genetic  variance  reported  by  allozymes  may  help  justify  increasing 
migration  among  MUs  through  Encroachment  Partnering. 

2.  Determine  the  degree  of  genetic  divergence  among  the  three  MUs.  MUs  are  defined 
genetically  as  independent  populations  in  which  the  rate  of  population  growth  is 
independent  of  rates  of  migration  (Moritz  1994;  Palsboll  et  al.  2007).  MUs  are  identified 
by  examining  the  level  of  genetic  divergence  among  local  populations  within  a 
landscape.  The  goal  is  to  identify  and  manage  units  in  a  way  that  will  help  maintain  the 
potential  to  adapt  to  environmental  conditions  while  not  necessarily  protecting  all  current 
adaptations  (Moritz  1994).  This  is  a  practical  objective  because  it  is  impossible  to 
identify  and  protect  the  genetic  basis  of  traits  that  will  be  adaptive  to  environmental 
conditions  in  the  future. 

Allozyme  data  will  likely  predict  only  one  population  unit  in  the  Onslow  Bight,  but  the 
microsatellite  data  may  predict  several  units,  because  their  spatial  distribution  will  be 
indicative  of  recent  changes  in  migration-drift  processes  due  to  their  high  mutation  rate 
and  selective  neutrality.  If  population  units  display  significant  divergence  in  allozyme 
markers,  then  genetic  differences  among  population  units  may  have  adaptive 
significance.  Therefore,  habitat  trading  that  increases  migration  would  still  increase 
population  viability  up  to  a  point,  but  high  rates  of  migration  may  cause  outbreeding 
depression.  While  we  doubt  this  will  be  the  case,  LEA  would  account  for  this 
contingency  by  including  three  genetic  services:  average  expected  heterozygosity  within 
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breeding  groups  (HS),  average  genetic  divergence  among  breeding  groups  within  an  MU 
(Dcs),  and  average  genetic  divergence  among  population  units  within  the  Onslow  Bight 
landscape  (Dst). 

3.  Generate  RCW  patterns  at  a  broader  scale  to  apply  Pattern-Oriented  Modeling. 

POM  could  be  used  to  calibrate  the  RCW  dispersal  submodel  derived  from 
radiotelemetry  studies  at  Fort  Bragg  (SI- 1471)  to  the  Onslow  Bight  population. 
Evaluation  of  parameter  uncertainty  in  the  RCW  SEPM  benefited  greatly  from  the  long¬ 
term  monitoring  program  on  MCBCL.  In  contrast,  the  only  consistent  pattern  across  all 
three  MUs  is  restricted  to  number  of  active  clusters  dating  back  to  1999.  Historical  data 
on  the  number  of  individuals  shared  among  the  three  MUs  is  virtually  absent,  except  for  a 
few  banded  birds  that  moved  between  Croatan  and  MCBCL.  Statistical  techniques  for 
summarizing  population  genetic  patterns  are  rapidly  developing  (e.g.,  Geneland:  Guillot 
et  al.  2005;  Structure:  Pritchard  et  al.  2000;  and  Bayesassnm:  Jehle  et  al.  2005).  We  will 
develop  techniques  to  link  genetic  patterns  observed  in  the  Onslow  Bight  to  those 
predicted  by  the  RCW  SEPM,  thereby  providing  an  empirical  basis  for  applying  POM 
across  the  Onslow  Bight. 

The  RCW  SEPM  with  the  empirically-derived  dispersal  submodel  will  be  used  to  predict 
baseline  levels  of  genetic  variance  within  and  among  population  units  expected  under  the  RCW 
recovery  plan.  However,  the  SEPM  can  also  be  used  to  estimate  the  historic  levels  of  genetic 
services  congruent  with  the  spatial  distribution  of  allozymes,  indicative  of  historic  rates  of 
migration  among  population  units.  Illustrating  how  Encroachment  Partnering  can  be  used  to 
move  rates  of  migration  closer  to  levels  observed  historically  will  provide  the  basis  for  shifting 
the  spatial  distribution  of  RCW  clusters  away  from  that  required  under  the  current  recovery  plan. 


4.2.2  Development  of  SEPMs  for  other  Species  (Bruggeman,  T  Wiegand  et  al.) 

It  has  recently  been  argued  that  habitat  fragmentation  is  a  poorly  defined  concept  (Fahrig  2003; 
Lindenmayer  and  Fischer  2007).  Researchers  have  not  used  a  consistent  approach  for 
quantifying  geometric,  or  structural,  changes  in  the  spatial  distribution  of  habitat.  Changes  in  the 
biological  system  are  also  quantified  in  different  ways  including  patch  occupancy,  total 
abundance,  migration  distances  or  rates,  and,  sometimes,  changes  in  spatial  apportionment  of 
genetic  variance.  The  field  of  landscape  ecology  has  also  been  hindered  by  summarizing 
patterns  observed  in  nature  and  then  proposing  theories  for  the  biological  processes  that 
generated  the  patterns  (i.e.,  process  is  inferred  from  pattern).  Inferring  process  from  pattern 
while  using  inconsistent  methodology  has  led  many  scientists  to  conclude  that  habitat 
fragmentation  effects  are  species-  and  landscape-specific,  which  can  not  be  generalized 
(Wiegand  et  al.  2005).  This  represents  a  large  technological  hurdle  for  military  installations 
wishing  to  shift  the  allocation  of  threatened  and  endangered  species  habitat  to  ensure  military 
readiness.  Below  we  present  some  surprising  results  that  will  provide  a  step  forward. 

Prior  to  beginning  research  on  our  SEED  grant,  the  Helmholtz  Centre  for  Environmental 
Research  -  UFZ  funded  collaboration  between  Douglas  Bruggeman  and  Thorsten  Wiegand  to 
establish  relationships  between  structural  and  biological  descriptors  of  habitat  fragmentation 
(Bruggeman  et  al.  submitted).  We  already  have  a  hypothesis  regarding  how  landscape  pattern 
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affects  biological  processes  for  RCWs  based  on  20  years  of  empirical  observations  on  military 
lands  (i.e.,  the  Letcher  model).  Below  we  describe  a  method  that  allowed  us  to  replace  that 
hypothesis  of  process  with  a  set  of  indices  describing  land  use  patterns  (i.e.,  pattern  was  inferred 
from  process).  We  argue  that  this  method  will  provide  a  generally  applicable  and  scientifically 
rigorous  approach  for  inferring  process  from  pattern  for  other  species  and  should  be  used  as  a 
starting  point  for  building  SEPMs. 

Methods:  Inferring  Pattern  from  Process 

The  RCW  SEPM  with  a  HCRW  dispersal  submodel  was  simulated  on  60  different  hypothetical 
landscapes  that  varied  from  clumped  to  fragmented.  Three  land  use  classes  were  employed 
including  old  growth  Longleaf  pine,  mixed  pine-hardwood,  and  non-forested.  The  HCRW 
submodel  assumed  that  dispersing  birds  avoid  crossing  non-forested  areas  and  showed  greatest 
preference  for  Longleaf  pine.  Both  the  first  and  second  order  effects  of  landscape  structure  were 
included  in  the  analysis  (Bailey  and  Gatrell  1995).  First  order  effects  summarize  the  mean 
density  of  RCW  habitat  across  the  entire  landscape.  We  used  three  densities  of  old  growth 
Longleaf  pine  or  RCW  breeding  habitat  (2.5%,  5%,  and  7%).  Second  order  effects,  which 
summarize  the  spatial  dependence  of  land  use  classes  at  small  spatial  scales,  were  summarized 
using  Ring  Indices  based  on  the  fundamental  the  physical  theory  of  point  pattern  phenomenon 
(Wiegand  et  al.  2005).  Ring  Indices  report  the  probabilities  of  finding  land  use  type  1  at  location 
X  conditioned  on  the  occurrence  of  finding  land  use  type  2  at  location  Y.  Where  r  =  X  -  Y,  or 
the  lag  distance  between  two  points,  Ring  Indices  are  symbolized  as  012(r).  In  our  study,  land 
use  type  1  represents  RCW  habitat,  2  represents  mix  pine-hardwood,  and  3  represents  non- 
forested.  For  example,  012(r)  represents  the  probability  of  finding  mixed  pine-hardwood  at  lag 
distance  r  from  RCW  habitat,  providing  an  estimate  of  the  contribution  mixed  pine-hardwood 
forest  makes  to  habitat  fragmentation  at  lag  distance  r.  In  contrast,  01  l(r)  provides  a  measure  of 
habitat  adjacency  at  small  spatial  lags  or  habitat  dispersion  at  large  spatial  lags. 

The  influence  of  first  and  second  order  effects  of  landscape  structure  on  biological  structure  was 
estimated  by  running  the  RCW  SEPM  for  100  years  to  achieve  equilibrium  between  landscape 
pattern  and  biological  processes.  Biological  results  were  summarized  as  N  (census  population 
size  or  number  of  adults  found  across  all  clusters),  F  (Wright’s  inbreeding  coefficient),  and  Fst 
(average  coancestry  within  a  breeding  group  due  to  population  subdivision). 

We  used  multiple  linear  regression  to  determine  if  the  dependent  variables  describing  population 
structure  could  be  predicted  by  independent  variables  describing  first  and  second  order  effects  of 
landscape  structure.  Analysis  began  with  null  regression  models  relating  population  structure  to 
known  first  order  effects.  For  abundance,  N  =  [3  *  fi  (Fahrig  2003),  where  1)  =  percent  of  RCW 
habitat.  For  genetic  indices,  G  =  P  *  log(N)  (Frankham  1996;  Wang  2004).  We  then  regressed 
the  residuals  from  these  null  models  against  Ring  Indices  01 1,  012,  013  for  lag  distances  of  1  to 
16  cells  (707  m  to  1 1.3  km).  Ring  Index-lag  distance  combinations  explaining  most  of  the 
residual  variance  was  assembled  into  a  “full  model”  and  AIC  model  selection  was  used  to  find 
the  most  parsimonious  multiple  regression  equation  (Burnham  and  Anderson  2002). 

Results 

Results  suggested  that  given  the  assumption  of  equilibrium  between  population  processes  and 
landscape  pattern,  the  RCW  SEPM  could  be  replaced  by  a  series  of  multiple  linear  regression 
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equations.  For  abundance  (Fig.  4. 1),  a  strong  relationship  with  first  order  effects  was  observed 
indicative  of  the  critical  role  of  breeding  habitat  for  predicting  population  size  (R2=0.98). 
However,  inclusion  of  Ring  Indices,  or  second  order  effects,  improved  model  fit  considerably 
based  on  reduction  in  AIC  scores  (AIC  Null  Model  =  -296,  AIC  Full  Model  =  -420). 

N  =  12.2*f  N  =  14.1  *fr  -  0.07*01 1(16)  -  0.06*012(1)  -  0.2*013(4) 

R2  =  0.98  R2  =  0.998 


Figure  4.1.  Predicted  census  population  size  (N)  using  the  RCW  SEPM  assuming  a  HCRW,  and 
N  predicted  by  the  null  and  full  linear  regression  models.  Diagonal  line  represents  the 
hypothetical  complete  agreement  between  the  two  modeling  approaches. 

The  strongest  fragmentation  effect  was  observed  due  to  the  presence  of  non-forested  areas  at  a 
distance  of  4  cells  (2.8  km)  from  breeding  habitat  (013(4)).  The  perceptual  distance  of  RCWs 
was  set  at  5  cells  for  these  simulations.  Therefore,  the  presence  of  non-forested  areas  at  the  edge 
of  a  fledgling’s,  helpers,  or  floaters  perceptual  distance  decreased  N.  Similarly,  the  presence  of 
mix  pine-hardwood  forest  adjacent  to  breeding  habitat  decreased  population  size  (012(1)).  The 
dispersion  of  breeding  habitat  across  the  landscape,  as  would  be  expected  when  01 1(16)  is  large, 
also  decreased  population  size.  The  influence  of  012(1)  and  011(16)  highlight  the  importance 
of  habitat  adjacency  included  in  the  model  structure.  The  model  rules  predict  that  more  helpers 
will  be  retained  under  a  clumped  distribution  of  habitat  and  that  these  helpers  increase 
recruitment  (Letcher  et  al.  1998;  Pasinelli  and  Walters  2002). 

Prediction  of  Wright’s  inbreeding  coefficient  (F)  was  more  difficult  (Fig.  4.2).  Total  population 
size  only  accounted  for  59%  of  the  variance  in  inbreeding.  However,  inclusion  of  the  ring 
indices  improved  model  fit  considerably,  increasing  R2  to  0.89  (AIC  Null  Model  =  -242,  AIC 
Full  Model  =  -321) 
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F  =  1.7  —  1.6*log(N)  F  =  0.98  -  1.0*log(N)  + 0.21*012(1) -0.3*013(4) 

R2  =  0.59  R2  =  0.89 


Figure  4.2.  Prediction  of  inbreeding  coefficient  (F)  by  the  RCW  SEPM  assuming  a  HCRW,  and 
by  the  null  and  full  linear  regression  models.  Note,  F  values  are  large  because  the  model  was  run 
for  100  years  without  the  addition  of  immigrants.  Diagonal  line  represents  the  hypothetical 
complete  agreement  between  the  two  modeling  approaches. 

Unlike  population  size,  a  considerable  portion  of  inbreeding  can  only  be  explained  when 
including  habitat  fragmentation  effects.  While,  011(16)  dropped  out  of  the  equation  during  AIC 
model  selection  procedures,  the  other  Ring  Index-lag  distance  combinations  are  consistent  across 
the  full  regression  models  for  total  population  size  and  inbreeding.  Therefore,  patterns  detected 
as  important  descriptors  of  biological  processes  hard- wired  into  the  SEPM  were  consistent  across 
these  statistical  summaries  of  population  structure. 

To  determine  if  the  Ring  Indices  might  detect  the  influence  of  landscape  friction  on  population 
structure,  we  also  ran  the  RCW  SEPM  assuming  straight  dispersal.  We  changed  the  model  rules 
to  decrease  floater  survival  based  on  the  number  of  steps  taken  in  non-forested  areas  to  simulate 
increased  predation  risks  hypothesized  for  small  forest  birds  (Belisle  and  Desrochers  2002). 
Comparing  results  across  the  HCRW  and  the  straight  dispersal  model  indicated  that  the  Ring 
Index-lag  distance  combinations  were  consistent  for  N  and  F.  However,  the  fitted  coefficients 
varied  in  a  manner  consistent  with  the  dispersal  rules  (data  not  shown).  Average  coancestry 
among  groups  (Fst)  is  a  more  sensitive  indicator  of  migration  rates.  Accordingly,  a  shift  in  Ring 
Index-lag  distance  combinations  was  observed  for  Fst-  Figure  4.3  indicates  that  there  is  no 
relationship  between  Fst  and  population  size  under  the  null  regression  model.  Inclusion  of 
second  order  effects  improved  model  fit  considerably! 
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Fst  =  0.82  -  0.09*log(N)  FST  =  0.65  +  0.33*log(N)  -  0.27*01 1(2)  -  0.09*012(8) 
R2  =  -0.003  R2  =  0.75 


0.5  0.6  0.7  0.8  0.9  1 

Predicted  by  RCWSEPM 

Figure  4.3.  Prediction  of  average  coancestry  due  to  population  subdivision  (Fst)  by  the  ROW 
SEPM  assuming  straight  dispersal,  and  by  the  null  and  full  linear  regression  models.  Note,  Fst 
values  are  large  because  the  model  was  run  for  100  years  without  the  addition  of  immigrants. 
Diagonal  line  represents  the  hypothetical  complete  agreement  between  the  two  modeling 
approaches. 

Population  geneticist  often  use  spatially-implicit  models  that  hold  second  order  effects  constant 
by  assuming  equal  rates  of  recruitment  within  each  population  unit  and  equal  rates  of  migration 
among  all  population  units.  Such  models  indicate  that  Fst  is  positively  related  with  log(N) 
(Wang  2004).  When  Ring  Indices  are  included  in  the  full  regression  model  this  theoretical 
expectation  holds.  This  suggests  that  variance  in  migration  rates  among  population  units  present 
in  a  spatially-explicit  model  are  accounted  for  by  the  Ring  Indices.  The  shift  in  Ring  Index-lag 
distance  combinations  is  also  supportive  of  this  conclusion.  Retention  of  012(8)  indicates  that 
when  mixed  pine-hardwood  forests  are  located  at  lag  distance  8,  Fst  will  be  lower.  This  result 
was  only  observed  when  floaters  moved  in  a  straight  line  and  suffered  a  reduction  in  survival 
when  crossing  non-forested  gaps.  Therefore,  landscapes  containing  higher  densities  of  forested 
areas  far  away  from  breeding  habitat  increased  migration  rates  among  groups,  thus  lowering  Fst- 
While  011(1  to  4)  accounted  for  a  significant  amount  of  the  residual  variance,  only  01 1(2)  was 
retained  after  AIC  model  selection.  Therefore,  01 1(2)  may  be  detecting  a  threshold  of  habitat 
adjacency  associated  with  a  decrease  in  Fst- 

Discussion  &  Future  Directions 

These  results  suggest  that  under  equilibrium  landscape  conditions,  an  SEPM  with  37  parameters 
(assuming  a  HCRW)  can  be  replaced  by  three  linear  regression  equations,  each  containing  at 
most  5  fitted  parameters.  We  stress  that  the  relationships  presented  are  very  likely  invalid  when 
landscape  structure  is  changing  over  time,  as  would  be  required  under  habitat  trading  programs. 
When  changing  landscape  structure  over  time,  SEPMs  are  is  still  preferred  because  ideally  they 
are  capable  of  capturing  interactions  between  landscape  pattern  and  biological  processes  that 
may  change  based  on  system  condition  in  the  previous  time  step. 


0.5  0.6  0.7  0.8  0.9  1 

Predicted  by  RCWSEPM 
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SEPMs  are  rule-based  models.  The  RCW  SEPM  incorporates  rules  that  relate  the  second  order 
spatial  structure  of  habitat  to  behavior  and  demography.  Integration  of  AIC  model  selection  with 
the  Ring  Index  approach  was  capable  of  detecting  spatial  thresholds  associated  with  these  rules. 
Are  other  species  that  limit  military  readiness  sensitive  to  second  order  effects  (i.e.,  habitat 
fragmentation)?  If  so,  Encroachment  Partnering  and  sustainability  of  biodiversity  would  benefit 
from  deriving  mechanistic,  rule-based  models  for  more  species  and  bases. 

We  propose  to  apply  this  approach  “in  reverse”  to  a  different  endangered  or  threatened  species. 

In  the  preliminary  study  just  presented  we  inferred  pattern  from  process,  thanks  to  the 
availability  of  an  SEPM.  By  perfonning  the  research  “in  reverse”  we  would  infer  process  from 
pattern  using  a  consistent  and  theoretically  rigorous  approach.  The  objective  would  be  to 
integrate  population  genetic  theory  with  the  theory  of  point  pattern  phenomenon  to  derive  rules 
for  a  new  SEPM.  To  do  so  would  require  conducting  a  landscape-scale  population  genetic  study 
as  outlined  for  RCWs  in  the  Onslow  Bight  above.  However,  the  studies  should  be  conducted  in 
landscapes  in  which  the  classification  of  land  cover  is  complete  or  nearly  complete.  Such  a 
study  would  benefit  from  identifying  a  species  that  that  has  limited  military  readiness  on  several 
bases.  We  have  not  identified  collaborators  or  bases  well  suited  to  the  approach  and  look 
forward  to  discussing  this  with  SERDP  in  the  near  future.  Our  strengths  lie  in  compiling 
information  in  a  manner  that  is  useful  for  determining  how  much  regulatory  relief  can  be 
achieved  by  trading  habitat  at  a  landscape-scale. 
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Appendix  A 


Apportioning  Genetic  Variance  over  Space 

Estimates  of  biodiversity  service  flows  in  the  form  of  genetic  variance  can  be  obtained 
using  a  method  derived  by  Nei  (1973).  Nei’s  method  estimates  the  frequency  of  heterozygous 
genotypes  based  on  the  frequency  of  alleles  observed  in  the  metapopulation.  This  method 
reports  the  heterozygosity  expected  if  the  metapopulation  were  in  Hardy-Weinberg  Equilibrium 
(i.e.,  alleles  are  assorting  in  offspring  assuming  a  random  mating  process).  This  approach  uses 
data  on  allele  frequencies  within  and  differences  among  subpopulations.  Nei  defines  genetic 
diversity  as  the  probability  that  any  two  alleles  chosen  at  random,  one  from  each  of  two 
individuals  are  independent.  By  assuming  that  individuals  are  mating  at  random,  this  probability 
is  equivalent  to  the  expected  heterozygosity.  The  two  alleles  can  be  sampled  from  within  an 
individual,  from  two  different  individuals  within  the  same  subpopulation,  or  from  two 
individuals  from  different  subpopulations  within  the  same  metapopulation.  The  expected 
heterozygosity  of  an  individual  in  subpopulation  i  (number  of  subpopulations,  S=l,  2,  3,  ...i)  can 
be  calculated  from  the  frequency  of  allele  m  (number  of  alleles,  k  =  1 ,  2,  3,  . .  .m)  in  the 
subpopulation  given  as  Xm: 

m 

H,  =1-2X-  (A.  1 ) 

k= 1 

H;  equals  the  average  proportion  of  heterozygotes  in  the  subpopulation.  Summarized  another 
way,  Hi  equals  the  average  genetic  diversity  of  an  individual  in  subpopulation  i.  If  k=2,  then  the 
maximum  value  of  Hi  is  0.5,  meaning  that  the  average  individual  is  heterozygous  at  that  loci. 
When  k>2,  H;  may  approach  unity  if  the  frequency  of  alleles  are  similar  to  each  other,  and  H; 
will  approach  zero  when  the  frequency  of  alleles  are  skewed.  The  probability  that  two  genes  are 
different  when  chosen  at  random  from  individuals  in  different  subpopulations  (e.g., 
subpopulation  i  and  j,  H,,)  is  calculated  as: 

m 

H,j=  1-EV-W.-  <A-2) 

k= 1 

The  average  expected  heterozygosity  within  each  subpopulation  (Hs)  can  be  estimated  by 
adding  a  weighting  factor  to  correct  for  differences  in  subpopulation  size  (w,  =  n,/N)  as: 

^  i  m  A 

ns=i-  •  <A3> 

V  5=1  k=\  J 

The  average  gene  diversity  among  subpopulations  (genetic  divergence,  Dst)  can  be  calculated  as 
the  average  differences  in  heterozygosities  between  all  pairwise  subpopulations: 


f  J 

Dst  =YjHw>w}  11 


Hi+Hj 


(A.  4) 


The  product  of  Hy  -  (H;  +  Hj)/2  is  Nei’s  minimum  genetic  distance  between  two  populations 
(Dm).  The  total  gene  diversity  of  the  metapopulation  can  now  be  estimated  by  the  separate 
contribution  of  intra-  and  inter-subpopulation  genetic  variances  (HT): 

H j.  =  H $  +  D .  (A. 5) 

Or  Hj  can  be  estimated  directly  from  allele  frequencies: 
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H  I  £  t,w,X„ 


(A.6) 


k= 1  V  *=1  / 

These  equations  do  not  account  for  the  statistical  bias  that  results  when  allele  frequencies  are 
estimated  from  a  subset  of  individuals  in  a  subpopulation.  Adjustments  to  these  equations  to 
account  for  sampling  bias  have  been  summarized  in  Nei  and  Chesser  (1983). 
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Camp  Lejeune.  Presentation  to  the  Onslow  Bight  Partnership.  April  28,  2006.  Habitat  Trading 
for  Red-Cockaded  Woodpeckers:  Incorporating  Landscape  Processes  and  Providing  Incentives 
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